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CHAPTER  I 


INTRODUCTION 


1.1  OBJECTIVES  AND  SCOPE 

The  objective  of  this  dissertation  is  to  develop  a  method  for 
determining  the  first  and  second  moments  of  the  dynamic  response  of 
structural  systems  with  random  parameters. 

The  randomness  in  the  response  of  structural  systems  is  due  to  two 
primary  factors.  The  first  of  these  is  the  randomness  in  the  loading 
which  has  been  discussed  extensively  in  the  past.  The  other  is  the  ran¬ 
domness  associated  with  the  structure  and  has  usually  been  ignored  in 
computing  the  moments  of  the  response. 

The  dissertation  will  develop  a  new  method  for  including  the  struc¬ 
tural  uncertainties  directly  in  the  computation  of  the  response  moments. 

The  systems  under  study  are  limited  to  those  with  linear,  elastic 
properties.  The  random  structural  properties  considered  are  the  stiff¬ 
ness  and  damping.  The  response  parameters  of  interest  are  the  deflec¬ 
tions  at  key  points  on  the  structure.  The  desired  results  are  the  mean 
and  covariances  of  the  deflections  as  functions  of  time. 

The  method  developed  will  not  require  knowledge  of  the  distribution 
of  the  structural  parameter  other  than  their  means  and  variances.  Sim¬ 
ilarly  the  distribution  of  the  response  will  not  be  computed.  For  the 
computation  of  first  passage  probabilities  the  response  will  be  assumed 
Gaussian. 

1.2  MOTIVATION 

The  randomness  in  structural  parameters,  such  as  stiffness, 
strength  and  damping,  has  long  been  realized  by  designers.  For  typical 
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structures  which  are  designed  to  resist  normal  fluctuations  in  environ¬ 
mental  loads,  such  uncertainties  may  be  safely  accounted  for  by  using 
conservative  values  for  the  stiffness  or  strength  rather  than  the  mean 
values.  However,  for  structures  which  must  resist  extreme  loads  such  as 
earthquakes,  accidental  blast,  or  weapons  effects,  such  a  conservative 
approach  is  often  economically  Impractical.  A  preferable  approach 
establishes  a  probability  of  failure  that  is  deemed  acceptable,  and  the 
structure  is  designed  to  insure  such  a  probability  is  not  exceeded. 

For  the  engineer  to  design  a  structure  with  a  specific  probability 
of  failure  three  things  are  required.  The  first  of  these  is  a  probab¬ 
ilistic  description  of  the  loading;  secondly,  the  probabilistic  charac¬ 
terization  of  the  structure  is  required.  Third,  and  the  most  difficult, 
is  a  method  for  computing  the  moments  of  the  response  and,  hence,  the 
probability  of  failure. 

The  probability  of  failure  is  the  probability  that  a  response  par¬ 
ameter  will  exceed  an  allowable  value  which  may  or  may  not  be  a  random 
variable.  When  the  allowable  value  is  deterministic  the  moments  of  the 
response  plus  the  assumption  or  computation  of  its  distribution  function 
will  be  sufficient  to  compute  the  probability  or  failure.  When  the 
allowable  value  is  a  random  variable  then  the  joint  distribution  of  the 
allowable  value  and  the  response  values  is  required. 

Typical  response  parameters  of  interest  in  computing  structural 
failure  include  strains,  stresses,  displacements,  velocities,  accelera¬ 
tions,  number  of  stress  reversals,  and  ductility  ratios.  Strains, 
stresses  and  ductility  ratios  are  usually  of  interest  when  examining  the 
possibility  of  overload  resulting  in  breaking  of  a  member.  Stress 
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reversals  may  be  of  interest  in  examining  fatigue  failure.  Displace¬ 
ment,  velocity  and  acceleration  are  of  interest  when  examining 
stability,  service  or  shock  isolation  type  failures.  For  example  when 
the  deflection  of  a  shock  isolated  piece  of  equipment  exceeds  available 
rattlespace  the  system  may  be  assumed  to  have  failed.  Similarly  a  tall 
building  may  suffer  a  service  failure  if  the  displacement  of  the  upper 
floors  in  high  winds  results  in  the  occupants  feeling  uncomfortable. 

This  type  of  failure  may  occur  when  the  stresses  and  strains  are  still 
well  below  a  dangerous  level. 

While  much  work  has  been  done  in  the  area  of  load  modeling  for 
earthquakes,  wind  and  blast  loads,  and  in  the  calculation  of  the  res¬ 
ponse  of  deterministic  systems  to  random  loads,  little  has  been  done  in 
establishing  techniques  for  computing  the  response  moments  in  the  pres¬ 
ence  of  randomness  in  the  structural  parameters. 

Recent  studies  described  in  the  literature  review  which  follows 
have  concentrated  on  quantifying  the  structural  uncertainties,  but  not 
on  their  effects  on  response  uncertainties. 

1.3  LITERATURE  REVIEW 

This  review  deals  with  studies  on  the  theory  and  calculation  of 
random  dynamic  structural  response.  Books  of  general  interest  will  be 
discussed  first  followed  by  papL~s  doling  with  the  quantification  of 
structural  uncertainties.  Papers  involving  the  computation  of  response 
characteristics  in  the  presence  of  structural  uncertainties  will  then  be 
discussed.  Finally,  some  papers  dealing  with  the  first  passage  and  peak 
values  of  structural  f-esponse  will  be  reviewed. 

The  theory  of  random  vibrations  of  linear  elastic  structures  with 
deterministic  mass,  damping  and  stiffness  properties  is  well  developed. 
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The  equation  of  motion  for  such  systems  when  they  are  mass  excited  is, 

I  [M]{X(t)l  +  [C]|X(t))  +  [K]{X(t))  =  {P(t)l  (1-1) 

where  [M]  =  System  Mass  Matrix 

I  [C]  =  System  Damping  Matrix 

[K]  =  System  Stiffness  Matrix 
{X}  =  System  Displacement  Vector 

j  {P(t)}  ■  System  Forcing  Vector 

When  {P(t)}  is  deterministic  {X(t)}  will  be  as  well, 

Lin  (1),  Crandall  (2),  and  Papoulis  (3)  all  address  the  stationary 
and  transient  response  of  such  linear  systems  for  both  stationary  and 
nonstationary  loadings,  Single-degree-of -freedom,  (SDOF),  systems  may 
ba  solved  using  the  system’s  impulse  response  function,  Multiple- 
^  degree-of -freedom,  (MDOF),  systems  are  analyzed  via  their  normal  modes. 

Expressions  for  the  autocorrelation  functions  and  first  passage  proba¬ 
bilities  of  the  response  ar»  developed  as  well  as  the  expressions  for 
^  the  mean  responses.  In  a  more  recent  work  Vanmar-ke  (4)  addresses  the 

extension  of  this  work  to  include  first  crossing  probabilities  for  both 
stationary  and  nonstationary  envelopes  about  the  process,  Elishakoff 
(5)  provides  an  excellent  text  for  quantifying  the  response  statistics 
for  bovh  static  and  dynamic  loads  as  well  as  a  discussion  of  the  Monte 
Carlo  method  of  analysis,  Soong  (6)  expresses  the  solution  of  differen 
tial  equations  with  random  coefficients  as  a  function  of  the  coeffi¬ 
cients.  The  mean  of  the  solution  at  any  time  may  then  be  obtained  by 


integrating  the  product  of  the  solution  and  the  joint  probability  den¬ 
sity  function  of  the  coefficients  over  the  range  of  the  coefficients. 

The  variance  may  be  obtained  in  a  similar  manner.  However,  the  result¬ 
ing  integrals  may  have  to  be  evaluated  numerically  at  each  time  step. 

Several  studies  have  addressed  the  quantification  of  structural 
uncertainties.  MacGregor,  Mirza  and  Ellingwood  (7)  compiled  information 
on  the  means  and  variances  for  concrete  strength,  steel  yield  strength, 
error  in  steel  placement  and  total  resistance.  Ellingwood  (8)  has  also 
examined  the  statistical  properties  of  reinforced  concrete  beam  col¬ 
umns.  These  studies  were  based  on  tests  of  existing  structures  and 
laboratory  specimens. 

An  analytical  study  was  performed  by  Grant,  Mirza  and  MacGregor  (9) 
to  examine  the  strength  of  concrete  columns  based  on  the  underlying  ran¬ 
dom  variables.  The  formula  for  the  strength  of  the  columns  was  evalua¬ 
ted  repeatedly  with  a  Monte  Carlo  procedure  to  determine  the  mean  and 
variance.  Ramsey,  Mirza  and  MacGregor  (10)  applied  a  similar  technique 
to  the  computation  of  deflections  of  concrete  beams.  Mlakar  (11)  used  a 
Taylor  series  expansion  about  the  mean  to  evaluate  the  mean  and  variance 
of  the  resistance  of  a  beam  column.  Thielen  (12)  used  a  similar  tech¬ 
nique  to  examine  the  behavior  of  concrete  beams  near  the  limit  state. 

The  consideration  of  structural  uncertainties  in  response  analysis 
IS  fairly  recent.  The  building  codes  have  recognized  strength  variation 
through  the  specification  of  strength  reduction  factors  and  working 
stresses  well  below  the  average  values.  The  field  of  reliability  and 
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safety  analysis  was  among  the  first  to  hoindle  the  uncertainties  explic¬ 
itly. 

Studies  within  this  field  have  concentrated  on  the  proper  design 
loads  and  strengths  to  use  for  structures  in  routine  circumstances.  Ang 
and  Cornell  (13)  investigated  the  influence  random  properties  such  as 
yield,  bond  and  shear  strength  have  on  structural  reliability  and  proper 
design  levels.  Ellingwood  and  Ang  (14)  proposed  procedures  for  risk 
based  design  as  did  Ravionda,  Lind  and  Siu  (15).  Moses  (16)  examined 
the  reliability  of  structural  systems. 

There  have  been  few  studies  concerned  with  the  calculation  cT  the 
moments  of  the  dynamic  response  of  a  system  with  random  stiffness  and 
damping.  Collins  and  Thompson  (17)  addressed  the  problem  of  random  eig¬ 
envalues  and  eignevectors  for  a  system  such  as  Equation  (1-1)  with  a 
random  stiffness  matrix.  Hasselman  and  Hart  (18)  developed  expressions 
for  the  modal  response  of  MDOF  systems  based  on  random  eigenvectors. 

In  both  of  these  papers  the  expressions  for  the  eigenvalue  means 
were  taken  to  be  the  same  values  that  would  resuH  from  the  determin¬ 
istic  expression  and  the  variances  were  evaluated  by  Taylor  series 
expansions  about  the  means.  Hisada  and  Nakgiri  (19)  used  a  Taylor 
series  expansion  about  the  means  to  quantify  the  effects  of  node  loca¬ 
tion  uncertainty  on  the  stiffness  matrix  and  resulting  stresses  and 
strains  for  the  static  case.  Handa  and  Anderson  (20)  developed  an 
approximation  for  solving  the  static  problem  by  separating  the  random 
variables  into  a  deterministic  variable  equal  to  the  mean  and  a  random 
variable  with  mean  zero  hut  a  variance  equal  to  the  variance  of  the 


original  random  variable.  However,  the  cross-correlation  between  the 
stiffness  and  deflection  is  ignored.  This  results  in  the  mean  response 
being  equal  to  the  response  of  the  system  when  the  system  parameters  are 
treated  as  being  deterministic  and  equal  to  their  means.  This  is  known 
to  be  incorrect. 

There  are  many  papers  related  to  first  passage  probabilities  and 
peak  values  for  stochastic  processes.  The  topic  is  also  well  covered  in 
many  of  the  basic  texts  listed  earlier. 

Ifi  a  classic  paper  Rice  (21)  developed  the  expressions  for  first 
passage  times  of  a  white  noise  process.  Karlin  and  Taylor  (22)  devel¬ 
oped  expressions  for  first  passage  probabilities  of  a  random  process 
represented  by  Brownian  motion  with  a  drift  in  the  mean.  Corotis,  Van- 
marke  and  Cornell  (23)  developed  expressions  for  the  first  passage  of 
nonstationary  random  processes.  Crandall,  Chandiramani  and  Cook  (24) 
used  a  nuinei^ical  simulation  technique  to  examine  first  passage  prob¬ 
abilities  of  the  response  of  a  SOOF  system.  The  system  was  linear, 
elastic  and  the  system  properties  were  deterministic.  The  loading  was 
repeatedly  applied  in  a  Monte  Carlo  style  and  the  statistics  of  the  res¬ 
ponse  including  the  expected  time  to  first  passage  beyond  a  barrier  com¬ 
puted.  A  number  of  types  of  barriers  were  investigated  as  were  the  zero 
and  stationary  start  conditions.  These  authors  also  developed  a  model 
where  the  sample  '^pace  is  divided  into  cells  and  the  probability  mass  is 
repeai'^'Jly  redistrlbcced  in  accordance  with  tneoretically  derived  tran¬ 
sition  probabilities,  krenk,  Madsen  and  Madsen  (25)  developed  a  method 
for  iTwdif\i'^g  stationary  crossing  frequencies  for  use  in  transient  res¬ 
ponse  envelopes.  Ang  (26)  developed  a  method  for  computing  the  response 
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moments  of  a  linear  system,  and  then  developed  first  passage  probabil 
ities  based  on  these.  Simple  nonlinear  elastic  structures  were 


examined. 


CHAPTER  2 


RESPONSE  OF  SDOF  SYSTEMS  WITH  RANDOM 
STIFFNESS  AND  DAMPING 


2.1  INTRODUCTION 

In  this  chapter  a  new  method  is  developed  for  analyzing  linear, 
elastic,  SDOF  systems  with  random  stiffness  and  damping.  The  system  may 
be  excited  by  either  a  deterministic  signal  or  random  process.  In 
either  case  the  response  of  the  system  will  be  a  stochastic  process. 

Such  processes  may  frequently  be  satisfactorily  characterized  by  their 
first  two  moments,  the  mean  and  variance.  For  the  processes  considered 
here  these  will  vary  as  functions  of  time. 

The  method  developed  involves  the  replacement  of  the  random  vari¬ 
ables  appearing  in  the  system's  equation  of  motion  with  two  new  vari¬ 
ables.  The  first  of  the  new  variables  is  deterministic  and  equal  to  the 
mean  of  the  original  random  variable.  The  second  variable  is  random 
with  mean  zero  and  has  a  variance  equal  to  that  of  the  original  random 
variable. 

This  substitution  has  been  used  by  Handa  and  Anderson  (20)  for 
analyzing  the  static  loading  of  structures  but  this  is  its  first  appli¬ 
cation  to  the  dynamic  response  of  structures.  Handa  and  Anderson  also 
neglect  the  correlation  between  the  randomness  in  the  system  parameters 
and  the  rando’  ness  in  the  response  introducing  an  error  in  the  mean 
which  they  acknowledge.  The  correlation  is  included  here  and  hence  the 
mean  generated  is  the  true  mean. 

The  results  of  the  method  are  the  mean  and  variance  of  the  system's 
displacement  and  velocity  as  functions  of  time.  The  correlation  between 
the  displacement  and  velocity  may  be  readily  computed  as  well. 


"T 
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The  details  of  the  transformation  of  the  random  variables  appearing 
in  the  equation  of  motion  are  explained  in  Section  2.2.  The  resulting 
transformed  equation  of  motion  is  solved  using  forward  differences,  the 
details  of  this  operation  are  explained  in  Section  2.3.  A  computer  code 
was  developed  to  perform  the  computations  and  its  basic  numerical  algor¬ 
ithm  is  provided  in  Section  2.4.  The  effects  of  some  higher  order  terms 
which  have  been  neglected  in  the  computation  of  the  variance  are  exam¬ 
ined  in  Section  2.5.  Several  numerical  examples  varifying  the  method 
are  provided  in  Section  2.6.  The  method  will  be  referred  to  as  the 
transformation  method. 

The  following  assumptions  are  made  concerning  the  system's  param¬ 
eters  and  response  characteristics, 

a)  The  mass  of  the  system  is  deterministic. 

b)  The  system  response  is  twice  mean  square  differentiable 
with  respect  to  time. 

c)  The  system's  damping  and  stiffness  are  random  variables 
whose  values  are  not  functions  of  time.  The  two  may  be 
correlated  or  uncorrelated. 

d)  When  the  forcing  function  is  stochastic  it  has  an  auto¬ 
correlation  function  with  one  factor  which  is  a  delta  function. 

e)  The  forcing  function  is  independent  of  the  system  parameters. 

2.2  TRANSFORMATION  OF  THE  EQUATIONS  OF  MOTION 

The  equation  of  motion  of  a  mass  excited  SDOF  system  is, 

MX  +  CX  +  KX  =  F(t)  (2-1) 

where,  M  =  System  Mass 

C  =  Damping  Resistance 
K  =  Stiffness 

X(t)  =  System  Displacement 
F(t)  =  Forcing  Function 
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When  K  and  C  are  deterministic  Equation  (2-1)  may  be  solved  exactly 
for  either  deterministic  or  stochastic  forciny  functions.  When  K  and  C 
are  random  variables  the  response,  X(t),  is  a  stochastic  process,  the 
mean  and  variance  of  which  may  be  evaluted  by  the  method  of  Soong  (6). 
However,  this  method  will  in  general  require  numerical  integration  at 
each  time  of  interest. 

To  develop  expressions  for  the  mean  and  mean  square  values  of  X(t) 
a  substitution  is  made  in  Equation  (2-1)  for  the  random  variables.  The 
substitution  replaces  the  original  variables  with  the  sum  of  two  other 
variables,  one  of  which  is  deterministic  and  equal  to  the  mean  of  the 
original  variable;  the  second  variable  is  random  with  mean  zero  and 
variance  equal  to  the  variances  of  the  original  variable.  The  substitu¬ 
tion  is  made  for  K,  C,  X(t),  X(t),  X{t)  and  F(t). 

The  deterministic  portion  of  the  response  X(t)  is  represented  by 
ii(t)  and  the  random  portion  by  6(t).  Hence  the  substitution  for  X{t) 
is, 

X(t)  =  u(t)  +  6(t)  (2-2) 

Taking  the  expected  value  of  Equation  (2-2)  yields, 

E(X(t))  =  E(u(t))  +  E(6(t))  (2-3a) 

=  u(t)  +  0  .  (2-3b) 

To  develop  an  expression  Tor  the  mean  square  value  of  the  response  one 
squares  Equation  (2-2)  which  yields. 


X(t)2  =  y(t)^  +  2u(t)6(t)+6(t)2  . 


(2-4) 


"■u 
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Taking  the  expected  value  of  Equation  (2-4)  yields  an  expression  for  the 
mean  square  response  of  X(t) 

E(X(T)2j  =  y(t)2  +  2ii(t).0  +  E(6(t)2)  (2-5a) 

E(X(t)2}  =  .,^t)2  +  E(a{t)2)  (2-5b) 

The  variance  of  X(t)  is  given  by, 

Var(X(t))  =  E{X{t)2j  -  E[X't))2  (2-6) 

Substituting  the  expressions  in  Equation  (2-bb)  for  E(X(t)^}  and 
Equation  (2-3b)  for  £{X(t))  yields, 

Var(X(t))  =  E(6(t)2)  (2-7) 

Similar  substitutions  may  be  developed  for  the  derivatives  of  X(t) 
by  differentiating  Equation  (2-2).  This  yields  for  the  velocity, 

X(t)  =  ;(t)  +  (5(t)  (2-8) 

And  for  the  acceleration, 

X(t)  =  p(t)  t  V(t)  (2-9) 

Taking  the  expected  value  of  Equation  (2-8)  yields  an  expression 
tor  the  mean  of  the  velocity, 

E(X(t))  =  E(u(t))  +  E(5(t))  (2-10) 

Since  u(t)  is  deterministic  Its  derivatives  are  as  well,  hence, 

£U(t))  =  u(t)  (2-11) 

Since  X(c)  and  thus  6{i)  have  been  assumed  mean  square  differentiable 


the  expectation  and  differentiation  operati'otis  may  be  interchanged.  And 
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since  E(6(t))  is  zero  for  all  times  the  expected  value  of  the  change  of 
6(t)  with  time  must  also  be  zero,  yielding, 

E(6(t))  =  0  .  (2-12) 

When  the  relation  in  Equations  (2-11)  and  (2-12)  are  substituted  into 
Equation  (2-10)  the  resulting  expression  for  the  mean  of  the  velocity 

■*5, 

E(X(t))  =  v{t)  .  (2-13) 

Squaring  Equation  (2-8)  and  taking  the  expected  valu'*  provides  the 
expression  for  the  mean  square  value  of  the  velocity, 

E(X(t)2)  =  u(t)2  +  E(5(t)^)  .  (2-14) 

Similar  expressions  for  the  mean  and  mean  square  values  of  the 
acceleration  may  be  developed  from  Equation  (2-9),  The  expected  value 
of  the  acceleration  is, 

E(X(t))  =u{t)  .  (2-15) 

The  mean  square  value  of  the  acceleration  is  given  by, 

E(X(t)2)  =  y(t)2  +  E(6(t)^]  *  (2-16) 

Representing  the  stiffness  K  as  the  sum  of  K[)  and  Kr  with  Kp 
deterministic  and  Kr  a  random  variable  with  mean  zero  yields, 

K  =  Kp  +  .  {2-1.’) 

faking  the  expected  value  of  Equation  (2-1/)  yields, 

E(K)  =  Kg  +  E(K.)  .  (;M3l 


14 


And  since  Kr  is  mean  zero, 

E(K)  =  Kp  .  (2-19) 

Squaring  Equation  (2-17)  to  develop  an  expression  for  the  square  of  K 
yields, 

=  kJ  +  2K^Kp+K^  .  (2-20) 

Taking  the  expected  value  of  Equation  (2-20)  yields 

E(k2)  =  +  E(k2)  .  (2-21) 

Upon  rearranging  one  has  the  expression  for  the  variance  of  K. 

Var(K)  =  E(k2)-E(K)^ 

=  E(K^)  . 

For  the  damping,  C,  let  Cq  and  Cr  be  the  deterministic  and 
random  portions,  respectively.  Then  following  the  development  used  for 
the  stiffness,  the  expressions  for  C,  the  mean  of  C,  and  the  variance  of 
C  are, 

c  =  Cjj  +  (2-23) 

E{C)  =  Cp  (2-24) 

Var(C)  =  E(c2)  .  (2-25) 

The  substitution  for  the  forcing  function,  F(t),  consists  of  the 
deterministic  portion  P{t)  and  a  random  portion  q(t).  Making  this  sub¬ 
stitution  and  following  Che  development  above  results  in  the  following 


(2-22a) 

(2-22b) 
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expressions  for  the  forcing  function,  F(t),  the  mean  forcing  function, 
P(t),  and  the  variance  of  the  forcing  function, 

F(t)  =  P{t)+q(t)  (2-26) 

E{F(t))  =  P(t)  (2-27) 

V.r(F(t))  =  E(q2(t))  .  (2-23) 

Substituting  the  results  of  Equations  (2-2),  (2-8),  (2-9),  (2-17), 
(2-23),  and  (2-261  into  the  equation  of  motion.  Equation  (2-1),  results 
in, 

M(u+6)  +  (Cp+C|^)(w+6)  +  (Kq+K|^)(u+6)  =  P(t)  +  q(t)  .  (2-29) 

lien  the  system  parameters  and  response  are  considered  correlated 
the  expected  value  of  Equation  (2-29)  is. 

Mu  +  CpU  +  Kpu  +  E(K^6)  +  E(C^5)  =  P(t)  .  (2-30) 

A  solution  for  u(t),  the  mean  response  of  the  system,  may  be  sought 
based  on  this  expression. 

To  develop  an  expression  for  5(t)  one  may  subtract  Equation  (2-30) 
from  Equation  (2-29)  yielding  the  result. 


M'6  +  Cq6  +  Kq6  +  JRT  +  KRT  =  q(t)  -  Kj^u  - 

(2-31a) 

where 

CRT  =  C^6-E(C|^6) 

(2-31b) 

KRT  =  Kj^6-E(Kj^6)  . 

(2-31C) 

The  term  CRT  represents  the  fluctuation  of  the  random  portion  of  the 
damping  force  about  its  mean.  The  fluctuation  of  the  random  portion  of 
the  stiffness  force  is  represented  by  KRT. 

When  the  random  portions  of  the  restoring  forces  are  considered 
smaT'  in  'mparison  to  the  deterministic  portions,  CpU  and  the 
fluctuations  may  be  considered  negligible  and  set  to  zero.  This  assump¬ 
tion  is  made  in  the  development  of  the  technique.  The  implications 
of  the  assumption  are  examined  in  Section  2.5  and  the  examples  in 
Section  2.6. 

When  the  CRT  and  KRT  terms  are  neglected  Equation  (2-31a)  becomes, 
m  +  Cp5  +  Kp6  =  q(t)  -  K^U  -  C^U  (2-32) 

The  response  characteristic  of  interest  in  this  analysis  is 
).  Furthermore,  the  equations  for  u(t)  and  6(t)  a,'e  coupled  in 
that  u(t)  appears  in  the  expression  for  6(t)  and  vice  versa.  The  sep¬ 
arate  the  equations  and  permit  the  computation  of  u(t)  and  E('S(t)^)  the 
forward  difference  technique  will  be  applied  as  explained  in  the  follow¬ 
ing  section. 

2.3  APPLICATION  OF  FINITE  OIFFERENCES 

Many  structural  mechanics  problems  involve  differential  equations. 
They  arise  in  both  static  and  dynamic  analysis  and  may  be  solved  dir¬ 
ectly  for  a  relatively  small  number  of  cases.  For  those  cases  which  ma/ 
not  be  solved  directly,  various  numerical  integration  schemes  have  been 
developed.  A  common  method  of  performing  the  integration  is  the  appli¬ 
cation  of  finite  differences.  This  technique  may  be  applied  to  a 
di tterentidl  equation  provided  the  function  and  its  derivatives  as  they 
appear  in  the  equation  are  continuous. 
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Consider  the  definition  of  a  derivative: 


df(y)  _  lim 
dy  '  Ay-i-o 


(2-33) 


When  Ay  approaches  a  small,  fixed  value,  s,  rather  than  zero  an 
approximation  to  the  derivative  of  f{y)  is: 


df(Y)  ^  f(y-^s)  -  f(y) 

W~  s 


(2-34) 


This  expression  is  the  basis  for  the  forward  difference  technique. 
If  one  considers  a  time,  t,  to  equal  nAt,  where  n  is  a  positive  integer 
and  At  is  a  small,  positive  time  constant,  the  expressions  for  fCi)  and 
its  derivatives  are. 


f(T)  =  f 


f(T)  = 


'n.r'n 


'n.2  -  ^ 


(2-35) 


(2-36) 


(2-37) 


For  further  information  on  the  finite  difference  approach  and  its  appli¬ 
cations  Wang  (31)  provides  a  good  introductory  discussion.  Rathe  (32) 
addresses  numerical  integration  schemes  in  consi de'^able  detail  as  well. 
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To  develop  the  relations  required  to  evaluate  y(t),  Equations 
(2-35),  (2-36),  and  (2-37)  are  applied  to  Equation  (2-30)  in  place  of 
u(t)  and  its  derivatives  resulting  in, 

^  'V2-2>'n+l''n)  (Vr"n'  S'“n>  “ 

At 


E(C,6)n  . 

(2-38) 

Combining  terms  and  solving 

for 

u„+2  yields. 

(P  -  E(Ko6) 
n  '  R  n 

(2-39a) 

where 

M  = 

=  M/At^ 

(2-39b) 

=  Cjj/At 

(2-39c) 

=  2M.Co 

(2-39d) 

Aj 

(2-39e) 

Equation 

(2-39a)  provides  i 

an  expression  by  which  u(t)  may 

be  calcu- 

lated  in  a  step-by-step  manner,  based  on  information  calculated  at  pre¬ 
vious  time  steps,  provided  that  information  can  be  developed.  The  main 
difficulty  lies  in  the  calculation  of  E(Kn6)  and  E(Cn6)  .  These  are 
the  terms  which  have  been  ignored  or  assumed  zero  in  many  previous  stud¬ 
ies.  Their  inclusion  here  represents  an  improvement  over  previous 
approxihk'tions.  The  means  for  calculating  these  terms  will  be  discussed 
shortly. 
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Applying  the  finite  difference  relations  of  Equations  (2-35), 
(2-36),  and  (2-37)  to  6(t)  and  its  derivatives  in  Equation  (2-32) 
yields. 


Rearranging  and  solving  for  6^+2, 


1 


=  S  '“In-  Vn-  "  Vnn"  ¥n>  '  '2-«> 


M 


Since  an  expression  for  'S(t)^  is  sought.  Equation  (2-41)  is  squared 
resulting  in, 

(‘'n  ■  ^‘’nVn  '  *  ^‘'n¥n+l 

M 

*  2<I„''2®ii  *  "Ri-n  "  2K|,C|,u^„ 


■  ^“n'^R^l^n+l  ■  ^“n'^R*2*n  *  ''R^ri 


-  ^'R%''2*n 


*  *1  «n4  *  *  *2  <  ) 


2  ,2 


(2-42) 


When  expected  values  are  taken  in  the  above  expression  some  terms 


on  the  right  side  will  be  zero.  Due  to  the  assumption  outlined  at  the 
end  of  Section  2.1,  terms  involving  a  product  of  Kr  and  qn  or  Cr  and 
qn  will  be  zero,  for  the  randomness  in  the  loading  is  independent  of  the 
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system  parameters.  Hence,  when  the  expected  value  of  Equation  (2-42)  is 
taken  the  result  is, 

Var(K)  VarC 
M 

.  Cov(K,C)  -  2u„AiE(K,6„,i) 

-  ^F,''2E(Vn>  - 

-2p„A2E(C,«„)  .  A^,E(6„4) 

+  A|E(5^„)  +  2AiA2E(«„,j5„)  +  2AjE(q„t„,l) 

Note  in  Equation  (2-43)  that  all  terms  which  involve  the  randomness 
of  the  system  parameters  are  multiplied  by  the  mean  response  or  its 
derivatives.  Hence,  when  the  system  has  a  zero  mean  response  the  vari¬ 
ance  and  mean  square  values  are  functions  of  the  randomness  of  the  exci¬ 
tation  only.  This  is  due  to  the  assumptions  regarding  the  terms  desig¬ 
nated  CRT  and  KRT  in  Equations  (2-31a),  (2-31b),  and  (2-31c).  These 
terms  are  discussed  .urther  in  Section  2.5. 

The  terms  Var(K),  Var(C),  Cov(K,C)  and  E(q^)  must  be  specified 
prior  to  commencing  the  solution.  The  terms  involving  E(6n+|‘^)  and 
B('Sn^)  0^  right  hand  side  are  the  results  of  Equation  (2-43)  at 
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previous  time  steps.  The  formulation  of  the  remaining  terms,  E(C„6  ), 
E(Kj^6^),  proceeds  as  follows.  Consider  Equation 

(2-41)  with  the  index  shifted  by  one.  This  provides  an  expression  for 
^n+1- 

M 

Multiplying  this  expression  by  Kr  yields 

Vn+1  "  ^'^R^n-1  “  "^R^n-l  ‘  '^R^R%-1 

*  W.  *  '^Vn-ll  •  <^-''5) 

Taking  the  expected  value  of  both  sides  results  in, 

E(Vnrt>  =  Y  (-Var(K)„„.j  -  p„.jCov(K.C) 

+  A^E(Kp6^)  +  A2E{y^_^))  .  (2-46) 

This  expression  may  then  be  readily  evaluated  in  a  step-by-step  manner 
from  the  zero  start  condition.  The  expression  for  E(Cp6^^^)  is  devel¬ 
oped  similarly  and  results  in, 

(-w„.iCov(K,C)  -l„.,Var(C) 
a  AjE{C„6„)  *  (2-47) 

The  expression  for  developed  in  a  similar  manner.  It 

is, 

r1 

•  A,E(S^„)  t  AjE(a„6„.,l) 


(2-48) 
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The  expression  for  E{q^6^^^)  is 


=i(E(q„q„.i)  -E(q„K,)u„.i 


(2-49) 


However,  since  the  loading  is  independent  of  the  system  parameters  and 
uncorrelated  with  itself  for  any  time  shift,  all  terms  in  the  above 
expression  are  zero.  Hence, 


Efq  6  =  E(q  6  )  =  0 

•^n  n+1  ^n  n 


(2-50) 


Since  these  terms  are  zero  the  randomness  in  the  response  due  to  the 
loading  and  that  due  to  the  system  parameters  are  independent.  Thus 
they  may  be  computed  separately  and  combined  merely  by  adding  the  two 
independent  variances. 

The  expressions  for  the  variance  of  the  velocity  may  be  developed 
from  Equations  (2-48),  (2-43)  and  the  finite  different  approximation  to 
the  velocity. 


6  ,  -  6 
n+1  n 


(2-37) 


Squaring  and  taking  the  expected  value  yields. 


,  E(6  J.)  -  2E(6  ^,6  )  +  E(6^’ 

,2  V  n+r  '  n+1  n'  '  n' 


(2-51) 


Hence  the  method  of  transforming  the  random  variables  by  a  shift 
a'>out  the  mean  provides  expressions  which  enable  the  analyst  to  evaluate 
the  means,  variances  and  covariance  of  displacement  and  velocity.  Note 
these  expressions  are  independent  of  the  distribution  of  the  system 
parameters  and  the  loading  function.  This  is  due  to  the  linearization 
of  the  solution  to  the  equation  of  motion  at  each  time  step. 

An  algorithm  for  computing  the  mean  and  variance  of  the  displace¬ 
ment  response  of  SDOF  systems  is  outlined  in  the  following  section. 

2.4  THE  COMPUTATION  ALGORITHM 

The  following  algorithm  was  developed  to  compute  the  mean  and  vari¬ 
ance  of  the  displacement  response  of  an  SDOF  system.  The  system  stiff¬ 
ness  and  damping  are  considered  random.  They  are  described  in  terms  of 
their  means  and  variances,  which  are  input  at  the  beginning  of  the  calc¬ 
ulation  as  is  their  correlation  coefficient.  The  randomness  in  the 
loading  is  described  in  terms  of  its  mean  and  variance  which  are  defined 
at  the  beginning  of  the  calculation  as  well.  These  may  vary  with  time 
provided  the  function  describing  the  variation  is  defined. 

The  algorithm  is  based  on  a  zero  start  condition  and  proceeds  as 
follows , 

1.  Read  input  parameters; 

K,  Var(K),  C,  Var(C),  Cov(K,C) 

M,  At,  and  problem  duration. 

2.  Initialize  all  variables  for  zero 
start  condition. 

3.  Set,  NT,  the  total  number  of  time  steps 
equal  to  the  problem  duration/At. 
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4.  Do  for  each  step  N=1,NT 

a)  Compute  the  mean  and  variance  of  the  loading 
function  based  on  predetermined  functions. 

b)  Compute  from  Equation  (2-48). 

c)  Compute  Equation  {2-39a). 

d)  Compute  E(5^^2)  Equation  (2-43). 

e)  Compute  (2-46). 

f)  compute  E(C|^'S^^2)  from  Equation  (2-47). 

g)  Write  desired  results. 

h)  Time  =  Time  +  At. 

i)  Go  to  step  4a  until 
Time  =  Problem  duration. 

This  computational  algorithm  is  implemented  in  the  computer  program 
RSDOF,  a  listing  of  which  is  included  in  Appendix  A.  Numerical  examples 
in  Section  2.6  use  the  code  to  varify  the  adequacy  of  the  method  as  well 
as  to  investigate  the  influence  the  random  system  parameters  have  on  the 
response  moments. 

2.5  THE  EFFECT  OF  HIGHER  ORDER  TERMS 

This  section  discusses  the  effects  of  the  assumptions  rade  in 
Section  2.2  regarding  the  fluctuations  in  the  random  portion  of  the 
system's  stiffness  and  damping  restoring  forces. 

Equation  (2-31a)  may  be  rewritten  as, 

M6  +  Cp5  +  Kp6  =  q(t)  -  -  Cj^y  -  CRT  -  KPT  f2-52) 

where  CRT  =  C[^6  -  E(C[^6)  (2-31b) 

KRT  =  -  E(Kj^6) 


(2-31c) 


25 


In  Section  2.2  KRT  and  CRT  were  assumed  to  be  zero.  Hence,  tney  do 
not  appear  in  the  development  of  the  finite  different  expressions  for 
Equation  (2-41)  and  E( ^2^’  Equation  (2-43)  in  Section  2.3.  In 
the  following  development  no  assumptions  regarding  the  values  of  CRT  and 
KRT  are  made  and  they  are  included  in  the  development  of  new  expressions 

for 

When  the  difference  expressions  are  substituted  for  6(v)  and  its 
derivatives  on  the  left  hand  side  of  Equation  (2-52)  the  result  is.- 

=  %  •  ^R^n-  '^R“n  - 


Solving  Equation  (2-53)  for  6n+2  yields. 


V2  =  ;  ■  ^R“n  -  *  W 

M 


(2-i)4) 


Note  the  difference  expressions  could  have  been  applied  to  the  6(t) 
terms  in  CRT  and  KRT  as  well,  but  this  would  only  complicate  the  algebra 
without  changing  the  results. 
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To  develop  an  expression  for  the  variance  of  the  response  one  first 
squares  Equation  (2-54)  resulting  in, 

V2  "  ^  ^‘^n  ■  ^^n'^R^n  "  ^^n^R^n  ^^n^'^n+1 

“  "  ^'■'z'S^'^n^n  S?’^n  "  ^^lSV]’‘n 

■  ^''2^^n^i  ^  ^l^n+1  ^^1^2'^n+l^n  ^2^n 

+  {CRT 2  +  2CRTKRT  +  KRT^  -  2q^CRT 

+  2u/rCRT  +  2u^C^CRT  -  2A^5^^^CRT 

-  2A26^CRT  -  2q^KRT  +  2Kj^  y^KRT 

+  2u^C^KRT  -  2Aj6^^jKRT  -  2A26^KRT})  (2-65) 

Equation  (2-55)  is  identical  to  Equation  (2-42)  developed  in 
Section  2.3  with  the  exception  of  the  terms  involving  CRT  and  KRT 
contained  within  the  brackets,  {  },  above.  Hence  when  one  takes  the 
expected  value  of  Equation  (2-55)  to  determine  the  variance  of  the 
response  the  result  is  Equation  (2-43)  plus  the  expected  value  of  the 
terms  within  the  brackets.  Thus  the  expected  value  of  the  terms  within 
the  brackets  provides  an  estimate  of  the  error  in  the  variance  of  the 
response  due  to  assuming  CRT  and  KRT  are  zero.  Designating  this  term 
REM  one  finds. 
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1 

REM  =  E  {CRT^  +  2CRT„KRT„  +  KRTJ 
n  t!2  *^11  n  n  n 

M 

*  ‘i^^n  -  %  -  'S'n+l  '  Vn> 

"  '<Vn  *  S“n  -  %  -  *l*nn  '  Vn>'^l  <''«> 

Evaluating  the  right  hand  side  of  Equation  (2-56)  term  by  term  yields, 

E(CRT2)  =  E[(C^6  -  E(C^6))2] 

=  Var(C^6) 

Likewise, 

E(KRT^)  =  Var(K,^6)  (2-58) 


And, 


E(KRTCRT)  =  E[(Cp6  -  E(C^6 ) )( K^6-E(K^6 ) ) ] 

=  Cov(K^6,C^6)  (2-59) 


For  the  remaining  terms  on  the  right  .hand  side  of  Equation  (2-56) 
involving  CRT  one  gets, 

S“n-  ■’n  '  ¥n*l  '  = 

V.  Vn»l-Vn>E‘V>nl  <2-60) 


Since  E(K|^),  E(C[^),  and  E(6^)  are  all  zero,  and  the  E(C[^<5)p  is 

deterministic  the  result  is, 

‘^R‘'n  ■  "n  '  ''iVl  '  ' 

=  »„E(VR«n)  *  in^C^R^n*  ' 

-  ''lE(6„,ii„C,)  -  A2E(6„i„C,)  (2-61 

The  corresponding  expression  for  the  terms  involving  KRT  in  Equation 
(2-56)  is, 

E((V„  t  C|,i„  -  q„  -  Aj5„^^  -  A25„)KRT„)  = 

'  *  “nE^VR^n*  ' 

-  ''l^'^nn'n^R*  -  ''2««n‘n^R) 

Combining  Equations  (2-57),  (2-53),  (2-59),  (2-61)  and  (2-62)  the 
resulting  expression  for  the  remainder  terms,  REM  is, 

REM  =  ^  (Var(K^5)^  +  ^  2Cov(K^6,C^6) J 

*  2lR„tE(Wn>  *  *  i„(E(K„C,6„)  .  E(C,^i„)] 

-  E'Wn'  -  E(AnVn)  -  ME(ViVr) 

*  E(6„S-„C,)])  (2-63 

This  expression  may  be  further  reduced  since  is  independent  of 
the  system  paraiiieters  Cj;^  and  Kj^  and  is  independent  of  and  6^  as 


I  ^  V  ^  I  \  K'  T  f  j  ■  j  w 


wwp^ipppi 
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Hence, 


E'^Vn’  =  E(q„)E{V„)  (2-64) 


And, 


''%V'  = 

However,  since  qp  is  mean  zero, 

'(%Vn)  =  ° 


(2-65) 


(2-66) 


And, 


E(q„v;.  =  0  (2-67) 

Thus  Equation  (2-63)  is  reduced  to, 

REM  -  (''ar(K^«)„  +  *  2Cov(K|,6,C|,5)^ 

H 

*  ^  EC^'n)! 

*  %lE<SVn!  ♦  E":r«„)1 

-  A2[E(«'Kr)  +  E(6„«„Cr)1  1  '  (2-68) 

Examination  of  the  tpriDs  on  the  right  hand  side,  of  Equation  (2-68) 
reveals  that  these  terms  provide  a  me^ns  for  the  uncertainties'  in  the 
system  parameters  entering  the  computation  of  the  response's  variance 


when  the  response  is  mean  zero. 


When  the  variation  of  the  response  is  small  compared  to  its  mean 
then  one  would  expect  that  the  value  of  Var(KR6)  would  be  small  com¬ 
pared  to  that  of  E{K^)u^.  Similarly  one  would  expect  that  the  Var(C|^6) 
would  be  small  compared  to  t(C^)u^  and  Cov(Kj^6,Cj^6)|^  would  be  small  in 
comparison  to  E(K^C^ Hence  one  would  expect  that  these  terms  have 
little  influence  on  the  response.  Similar  arguments  may  be  made  for 
discounting  the  rest  of  the  terms  in  Equation  (2-68)  as  they  all  have 
corresponding  terms  which  one  would  expect  to  be  much  larger  in  Equation 
(2-43). 

Several  of  the  following  examples  verify  this  assumption  by  compar¬ 
ing  the  results  of  the  computer  code  RSDOF  to  analytical  methods  which 
do  not  ignore  the  terms  included  in  Equation  (2-68). 

2.6  NUMERICAL  EXAMPLES 

This  section  includes  examples  of  the  transformation  technique 
developed  earlier  in  this  chapter  applied  to  the  solution  of  specific 
problems.  The  following  problems  are  covered. 

1)  Response  of  deterministic  random  systems  to  a 
white  noise  loading. 

2)  Response  of  a  random  system  to  a  deterministic 
constant  loading. 

3)  Response  of  a  random  system  to  a  deternr ni Stic  blast  loading. 

Tne  first  two  problems  were  used  as  a  tnetnod  of  checking  the 
results  of  L;'.e  technique.  These  examples  were  also  used  to  quantify  the 
effects  the  random  system  parameters  have  on  the  response.  The  third 
demonstrates  the  use  of  the  technique  in  a  practical  application. 
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EXAMPLE  1  -  White  Noise  Excitation  of  An  SDOF  System 

Consider  the  system  diagramed  in  Figure  2-1.  All  system  parameters 
are  unitless.  When  F(t)  is  white  noise  with  spectral  density  ij)  the  mean 
square  response  of  the  system  with  deterministic  parameters,  is  given  by 
Lin  (1)  as, 


E(X(t)2)  =  -Si-- 

2eu)Q  m^ 


exp(-2eu  t ) 

- H 

“D 


+  2(eojQSina)|jt)^  +  eu^si  n2iupt  ]  } 


(2-69) 


where 


2^ 


(2-70) 


(2-71) 


(2-72) 


As  the  response  reaches  its  stationary  value  Equation  (2-69) 


reduces  to. 


£s(X(t)2) 

2eu)^m^ 

0 


(2-73) 


To  model  the  response  of  the  system  one  first  has  to  develop  a 
suitable  model  for  the  white  noise  excitation.  A  wiiite  noise  random 
process  with  spectral  density  (j*  has  an  autocorrelation  function,  Ppp(i) 


Rpp  (t  )  --  2ti(|>S{  t) 


(2-74) 


where  6|rl  is  the  di rac  delta  function. 
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When  one  wishes  to  excite  a  structure  using  a  white  noise  input  and 
perforin  response  computations  in  a  discrete  time  framework  he  can  model 
the  input  as  a  frequency  of  independent,  identically  distributed  random 
variables  applied  to  the  structure  at  At  time  intervals.  This  is  an  IID 
(independent  identically  distributed)  random  process  (sometimes  referred 
to  as  a  band-limited  white  noise  random  process).  Let  {F.  =  F(jAt), 
j  =  0,  1,  2,  ...}  be  an  IID  random  process.  The  random  process  has  mean 
zero  and  autocorrelation  function. 


Rcr((j-k)at)  =  {  f  ^  (2-75) 

=  (2-76) 

where  5^  is  the  Kronecker  delta  (equal  to  zero  when  the  subscript  is 
nonzero,  and  equal  to  one  when  the  subscript  is  zero).  The  spectral 
density  of  this  random  process  is. 


S  (oj) 
FF 


Atq^ 

2ir 

0 


rt 


elsewhere 


If 

At 


(2-77) 


In  view  of  this,  the  IID  random  process  is  equivalent  to  the  white  noise 
random  process  if  is  chosen  as  follows, 


0 


2 


2tt 
"  At 


(2-78) 


Th’s  aefinition  for  insures  equivalence  between  the  white  noise 
and  IID  random  process  in  terms  of  their  spectral  densities.  However, 
equivalence  in  terms  of  the  structural  response  these  inputs  excite  is 
only  guaranteed  for  certain  values  of  At. 
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Specifically,  the  cutoff  frequency  of  the  I  ID  random  process,  ir/At, 
must  be  much  greater  than  the  natural  frequency  of  the  structure  being 
excited  in  order  for  the  IID  random  process  to  appear  similar  to  the 
white  noise  input,  as  far  as  structural  response  is  concerned.  This 
establishes  the  requirement, 

01 

At  «  —  (2-79) 

n 

where  oip  is  the  natural  frequency  of  the  structure  being  excited. 

To  address  the  issue  of  how  small  a  time  step.  At,  is  required  to 
insure  an  adequate  approximation  of  the  structural  response  character¬ 
istics  a  series  of  calculations  with  decreasing  values  c"  At  were  per¬ 
formed.  The  ratio  of  the  natural  period  to  the  time  step  is  termed  N. 
The  computations  were  run  for  the  case  of  C  =  20,  (e  =  0.2),  until  the 
response  became  stationary.  With  =  100/At,  K  *  500,  and  M  =  5,  the 
exact  solution  for  the  mean  square  value  of  the  response,  Es[X(t)^), 
as  defined  in  Equation  (2-73)  is  0.005.  The  values  for  Es(X(t)^)  com¬ 
puted  for  successive  values  of  N  are  plotted  vs  ln(N)  in  Figure  2-2. 

As  At  becomes  small  (that  is  as  N  becomes  large),  the  values 
obtained  from  the  computations  very  closely  approximate  the  exact  solu¬ 
tion,  being  in  error  by  less  than  1%  for  values  of  N  greater*  than  1500. 

To  examine  the  variation  of  the  error  with  time  the  results  of  com¬ 
putations  for  different  values  of  the  damping  parameter  C,  were  compared 
to  the  exact  solution  as  given  by  Equation  (2-69).  the  values  of  C 
chosen  for  the  comparison  correspond  to  values  of  e,  the  critical  damp¬ 
ing  ratio,  of  0.05,  0.10,  and  0.2.  The  computed  values  of  E(x(t)^)  are 
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plotted  versus  time  in  Figure  2-3.  The  relative  error,  E,  was  computed 
<  *  .time  step  by, 


E  - 


E(X(t)  ](-omp  ~  ^(X(t)  )gxac 


X  100 


(2-80) 


where  E(x(t)^)^Q^p  is  the  value  of  the  mean  square  response  as  computed 
and  ^(x(t)^]g^g^  is  the  value  for  the  mean  square  response  given  by 
Equation  (2-69).  The  error  was  found  to  be  monotonical ly  increasing 
with  time  tor  all  three  cases,  reaching  its  maximum  value  when  the  sys¬ 
tem  reached  the  point  of  stationary  response. 

The  values  for  E„,  ,  the  maximum  error  for  the  three  cases  are  lis- 
max 

ted  in  Table  (2-1).  The  maximum  error  was  found  to  increase  with  de¬ 
creasing  values  of  e.  However,  since  the  mean  square  response  goes  to 
infinity  as  e  goes  to  zero  this  is  not  surprising.  A  smaller  value  of 
At  would  result  in  smaller  errors  for  the  very  lightly  damped  system. 

If  one  now  considers  the  system  diagramed  in  Figure  2-1  to  have 
both  random  stiffness  K  and  damping  C  the  values  obtained  from  Equation 
(2-42)  for  the  mean  square  response  will  be  the  same  as  for  when  the 
system  is  deterministic.  This  is  due  to  the  assumption  made  in  Section 
2.1  regarding  the  variations  of  the  random  components  of  the  restoring 
forces.  Specifically  it  was  assumed  that, 

E(V)  =  C^6 

and  (2-81) 

E(Kj^6)  = 


TABLE  2-1 


MAX  ERROR  IN  THE  COMPUTED  VALUE  OF 
EfX(t)2)  FOR  VARIOUS  VALUES  OF  € 


€ 

E(X(‘)^)comp 

Max  error 

0.05 

0.020472 

2.55% 

0.10 

0.010128 

1.28% 

0.20 

0.005032 

0.64% 

2500  for  all  cases 
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These  assumptions  have  no  effect  on  the  mean  response  of  the  system 
but  do  result  in  the  contribution  of  the  randomness  in  the  stiffness  and 
damping  to  the  randomness  of  the  response  being  multiplied  by  a  factor 
containing  the  mean  response.  Hence  when  the  mean  response  is  zero  the 
effects  of  the  randomness  in  the  system  parameters  is  lost. 

To  answer  the  question  as  to  how  significant  this  error  is  one  may 
write  the  expression  for  the  mean  square  response  of  the  stationary  sys¬ 
tem  as  a  function  of  K  and  C.  Rearranging  Equation  (2-73)  results  in, 

Es(X(t)2)  =  (2-82) 

For  purposes  of  illustration  one  may  consider  K  random  and  C  deter¬ 
ministic,  resulting  in  Equation  (2-82)  being  a  function  of  the  random 
variable  K.  Taking  the  expected  value  of  Equation  (2-82)  results  in, 

E(E5(X(t)2))  =  E({)  (2-83) 

And  since  the  expectation  of  an  expected  value  is  the  expected  value  the 
left  hand  side  may  be  rewritten  and  Equation  (2-83)  becomes, 

Es(X(t;2)  =  ^  E(  -^  )  (2-84) 

The  E(^)  depends  on  the  distribution  of  K.  When  the  m(;an  stiff¬ 
ness  Kp,  is  equal  to  500  and  the  coefficient  of  variation  6|^  is  0.15, 


and  K  is  uniformly  distributed  with  upper  bound  and  lower  band 
E(i)  is  given  by. 


crii  =  fU  i  _J_ 

KJ  K  SPAN  “ 


12-85 


where  SPAN  =  Ky  -  Kl  and  is  given  by 


Var  K  =  i  {SPAN)2 


(2-86 


And  since, 


Var  K  =  (0.15) , 
Var  K  =  5625  ^  J 


(2-87 


The  span  of  the  distribution  is. 


SPAN  =  259.8 


(2-88 


Yielding  values  of  and  the  upper  and  lower  limits  of  the  distril 


tion  of, 


629.9  j 
=  370.1 


(2-89 


The  resulting  value  for  E(  -^  )  is  thus, 

E(  i  )  =  0.002046  =  ^ 


(2-90 


this  differs  from  —n—  by  only  2.35%.  Thus  the  true  value  for  the  sta- 
tion  any  mean  square  response  of  a  white  noise  excited  system  with 
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random  stiffness  is  underpredicted  by  only  about  2.35%  for  the  case 
illustrated.  A  suitable  rule  of  thumb  would  be 

^y)  (2-91) 

where  6^^  is  the  coefficient  of  variation  for  K  and  is  the  mean  value. 
This  example  provides  considerable  reassurance  that  the  assumptions  made 
in  Section  2.1  will  produce  only  a  small  ei .  >jr. 


■‘I 
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EXAMPLE  2  -  Comparison  of  the  Results  of  the  Transformation  Method  to 
the  Exact  Solution  for  the  Response  of  an  SDOF  System  with 
Random  Parameters 

The  objective  of  this  example  is  the  verification  of  the 
transformation  method  by  comparing  its  results  with  those  obtained  by 
numerically  integrating  an  exact  solution  for  the  mean  and  variance  of 
the  system  response. 

The  system  under  consideration  has  the  following  unitless  mean 
values, 

Kjj  =  500 
e  =  0.1 
M  =  5 

and  zero  initial  conditions. 

The  loading  function  used  for  the  comparison  was  deterministic. 

The  load  was  instantaneously  applied  at  time  0  and  remained  in  place 
indefinitely.  The  loading  function  is  diagramed  in  Figure  2-4. 

The  system  response  to  such  a  loading  is, 

c  eo) 

x(t)  =  (1  -  exp{-e!i)Qt)(— -  sinajjjt  +  cosu>pt))  (2-92) 

where  =JY  (2-93) 

u)p  =  (2-94) 

The  first  comparison  consists  of  evaluating  the  response  of  the 
SDOF  system  using  the  code  RSDOF  with  the  randomness  in  the  system  par¬ 
ameters  set  to  zero.  The  results  of  thii>  computation  were  plotted  with 
values  obtain-  from  Equation  (2-92)  in  Figure  2-5.  As  can  be  seen  in 
the  figure  t  plots  are  indistinguishable. 
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For  the  second  comparison  the  stiffness  was  treated  as  a  random 
variable  with  a  coefficient  of  variation  of  0.1.  To  generate  the  exact 
solution  for  the  mean  and  variance  of  X(t)  Equation  (2-92)  was  treated 
as  a  function  of  K  at  any  time  t.  Thus, 

X{t)  =  G{K) 


E(X(t))  =  E(G(K))  =  _^r  G(K)P^  dK 

(2-95) 

E(x2(t)]  =  E(G^(K))  g2(K)P^  dK 

(2-96) 

Var(X(t))  =  E(g2(k))  -  E(G(k))2 

(2-97) 

where  Px  is  the  probability  density  function  for  K. 

One  notes  that  the  exact  values  for  E(X(t))  and  Var[X(t))  are 
functions  of  the  probability  distribution  function  of  K.  This  informa¬ 
tion  Is  not  included  in  the  equations  of  the  transformation  method  since 
the  method  replaces  the  derivatives  with  linear  approximations.  The 
means  and  variances  of  the  resulting  linear  equations  are  independent  of 
the  distributions. 

The  integrations  of  Equations  (2-95)  and  (2-96)  were  performed  for 
two  different  distributions  of  the  random  variable  K.  The  first  of 
these  was  the  uniform  distribution.  The  probability  density  function 
for  K  is  then, 


P 


K 


K|^  <  K  < 


(2-98) 


0 


el spwhere 
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where  K^j  is  the  upper  limit  of  integration  and  the  lower.  Values  of 
and  were  413.4  and  586.6  respectively.  These  yield  a  mean  value 
of  K  of  500  and  a  coefficient  of  variation  of  .1. 

The  numerical  integration  was  performed  using  Newton-Quotes  algor¬ 
ithm  available  via  the  QUANC8  subroutine  available  on  the  University  of 
New  Mexico's  VSPC  computer  system. 

The  resulting  means  and  standard  deviations  for  approximately  2.5 
cycles  of  response  are  plotted  in  Figures  2-6  and  2-7  respectively.  The 
means  are  indistinguishable.  The  transformation  method  appears  to  ser¬ 
iously  undershoot  the  troughs  of  the  standard  deviation  curve.  However, 
since  the  response  at  or  near  the  peaks  is  of  primary  concern  to  the 
engineer  this  is  not  considered  a  serious  shortcoming. 

For  a  second  comparison  the  stiffness  K  was  assumed  to  have  a  trun¬ 
cated  normal  distribution.  The  probability  density  function  for  K  is 
then, 


P 


K 


1.00067 


(2-99) 


for  <  K  <  Ky 
=  0  elsewhere 


The  means  are  once  again  indistinguishable,  and  are  plotted  in 
Figure  2-8.  The  standard  deviations  are  plotted  in  Figure  2-9.  Once 
again  the  transformation  method  undershoots  the  troughs  of  tne  standard 
deviation  vs  time  curve. 

Based  on  the  abovf  information  the  transformation  method  provides 
an  adequate  approximation  for  the  mean  and  standard  deviation,  or  vari¬ 
ance,  of  an  SDOF  system  with  random  parameters. 
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EXAMPLE  3  -  Evaluation  of  the  Influence  of  Random  Stiffness  and  Damping 
on  the  Response  Statistics  for  an  SDOF  System 

This  example  consists  of  a  set  of  computations  of  the  mean  and 
standard  deviation  for  an  SDOF  system.  The  system's  mean  stiffness  and 
mass  were  held  constant.  The  loading  function  was  the  same  as  in 
Example  2  and  is  shown  in  Figure  2-4.  The  mean  damping  ratio,  the 
standard  deviations  for  the  stiffness  and  damping,  and  the  correlation 
coefficient  between  the  stiffness  and  damping  were  varied  in  individual 
cases.  The  values  used  for  these  parameters  in  each  case  are  shown  in 
Table  2-2. 

The  results  of  cases  1  and  2  and  the  deterministic  case  were  used 
to  examine  the  influence  of  the  randomness  in  the  stiffness  on  the  ran¬ 
domness  of  the  response.  The  resulting  values  for  the  mean,  u(t),  and 
the  standard  deviation,  o(t),  for  each  case  are  plotted  in  Figures  2-10, 
and  2-11.  The  deterministic  case  was  discussed  in  Example  2  and  its 
mean  response  is  plotted  in  Figure  2-5.  One  notes  that  the  value  for 
the  peak  response  increases  as  the  standard  deviation  of  the  stiffness 
increases.  Furthermore  the  standard  deviation  curve  lags  behind  the 
response  curve  though  they  are  of  similar  shape.  The  trend  of  increas¬ 
ing  peak  response  is  evident  by  comparing  the  values  for  the  peak  res¬ 
ponse  obtained  for  cases  1  and  2  and  the  deterministic  value.  However, 
the  variation  between  these  cases  is  smaller  than  it  would  be  for  the 
static  case,  hence  the  variation  in  the  stiffness  has  a  smaller  effect 
on  the  peak  dynamic  response  than  it  does  on  the  static.  Furthermore 
the  values  obtained  from  the  code  are  at  preselected  times  based  on  the 
period  of  the  response  as  determined  by  its  mean  parameters,  when  the 
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system  has  random  parameters  there  can  be  a  slight  phase  shift  resulting 
in  the  peak  values  printed  being  off  by  one  or  two  time  steps  from  the 
true  peaks. 

The  peak  response  for  each  of  the  cases,  a|^  =  0,  50,  100  are  listed 
in  Table  2-3. 

The  results  of  cases  4,  5,  and  6  allow  one  to  examine  the  influence 
of  the  randomness  in  the  damping  on  the  randomness  of  the  response.  The 
means  and  standard  deviations  for  these  cases  are  plotted  in  Figures 
2-12,  2-13,  and  2-14  respectively. 

By  comparing  the  results  of  cases  1  and  4  one  notes  that  the  values 
obtained  for  crx{t)  resulting  from  a  random  stiffness  are  an  order  of 
magnitude  larger  than  those  obtained  for  random  damping.  Cases  2  and  5 
show  the  same  effect  for  larger  variations  in  the  stiffness  and  damp¬ 
ing.  However,  one  must  note  that  the;>e  results  are  for  a  system  with 
damping  which  is  only  10%  of  the  critical  value.  When  the  damping  is 
increased  to  50%  of  critical  as  in  cases  3  and  6,  plotted  in  Figures 
2-15  and  2-14  respectively,  the  resulting  values  for  oylt)  resulting 
from  randomness  in  the  stiffness  and  damping  are  still  not  of  the  same 
magnitude. 

From  these  results  one  concludes  that  for  lightly  damped  systems 
which  do  not  have  a  zero  mean,  the  effects  on  response  of  randomness  in 
the  damping  are  minimal  compared  to  those  resulting  from  a  randomness  in 
the  stiffness.  For  the  case  of  white  noise  response  the  effects  of 
randomness  in  K  and  randomness  in  C  will  be  the  same  when  the 
coefficients  of  variation  of  K  and  C  are  equal. 
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Furthermore  while  the  randomness  in  the  stiffness  effects  the  res¬ 


ponse  for  its  entire  duration,  the  effects  of  the  random  damping  dissi¬ 
pate  as  the  system  converces  to  its  static  deflection  which  is  indepen¬ 
dent  of  the  damping.  However,  when  the  loading  is  transitory  the 
effects  of  the  random  stiffness  will  also  dissipate  as  the  system  re¬ 
turns  to  its  initial  state. 

The  results  cf  cases  7,  8,  and  9  allow  one  to  assess  the  impact  of 
correlation  between  the  stiffness  and  damping.  In  these  cases  the  cor¬ 
relation  coefficient  between  the  stiffness  and  damping,  Pkd,  was  given 
the  values  of  0,  1,  and  -1,  respectively. 

The  difference  between  these  cases  were  too  small  to  be  effectively 


plotted,  however,  case  8,  P^c  =  1.  produced  the  largest  values  of  both 


Uj^(t)  and  aj^(t).  Case  9,  P||,q 
and  ax(t). 


1,  produced  the  smallest  values  of  p{t) 


Bsed  on  these  results  one  concludes  that  ax(t)  at  time  of  peak 


response  increases  with  increasing  correlation  between  the  stiffness  and 


damping.  The  values  for  the  mean  and  standard  deviation  at  the  time  of 


peak  response  are  provided  in  Table  2-4  as  an  illustration. 
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EXAMPLE  4  -  Blast  Loading  of  a  Simply  Supported  Beam 

This  example  illustrates  the  application  of  the  methodology  devel¬ 
oped  in  the  previous  sections  to  a  problem  of  practical  interest. 

The  structure  considered  is  a  simply  supported  reinforced  concrete 
beam.  The  loading  is  a  simplified  approximation  to  a  blast  loading. 

The  structure  and  its  loading  time  history  are  shown  in  Figure  2-16. 

A  deterministic  solution  using  standard  analysis  procedures  is  pro¬ 
vided  in  Biggs  (33).  The  peak  deflection  resulting  from  the  determin¬ 
istic  and  stochastic  solutions  are  compared. 

When  modeling  sucn  a  structure  with  an  SDOF  system  the  stiffness 

is, 


Kp,  =  °  ^  -  (2-100) 

for  the  present  case  this  yields  a  value  of  Kq  of  2.08  x  10^  Ib/in. 

The  stiffness  is  assumed  to  have  a  coefficient  of  variation  of  0.15. 

This  results  in  the  variance  of  K,E(K|<^),  being  9.734  x  10  ^°  Ib^/in^. 

The  damping  was  assumed  deterministic  and  equal  to  8%  of  critical. 
The  results  of  the  stochastic  computation,  the  mean  of  the  midspan 
deflection  and  its  standard  deviation  are  plotted  in  Figure  2-17.  The 
maximum  value  obtained  for  the  mean  deflection  was  0.2450".  Biggs  esti¬ 
mated  an  upper  bound  of  0.2726"  for  the  undamped  structure. 

The  key  question  of  interest  for  the  engineer  is,  "what  is  the  pro¬ 
bability  of  the  peak  response  exceeding  a  desired  level,  b?".  This 
question  is  addressed  in  Example  2  of  Chapter  4  where  the  passage  prob¬ 
abilities  for  this  problem  are  computed. 


CHAPTER  3 


RESPONSE  OF  MDOF  SYSTEMS  WITH 
RANDOM  STIFFNESS  AND  DAMPING 


3.1  INTRODUCTION 

In  this  chapter  the  transformation  method  developed  in  Chapter  2  is 
expanded  to  MDOF  sytems.  This  introduces  a  number  of  complications  in¬ 
volving  the  correlation  between  the  random  variables  describing  the 
system  and  the  response  parameters.  Expressions  for  the  general  case 
where  the  correlation  between  system  parameters  may  be  user  specified 
and  for  a  simpler  case  of  individual  terms  of  the  system  matrices  being 
uncorrelated  are  developed.  The  latter  case  is  used  in  the  comptuer 
code  developed  simply  to  reduce  computation  times. 

The  assumptions  regarding  the  stochastic  nature  of  the  displacement 
response,  X(t),  and  of  the  loading,  F(t),  for  the  SDOF  system  of  Chapter 
2  to  apply  to  {X(t)},  the  vector  displacement  response,  and  {F(t)},  the 
forcing  function  vector  as  well. 

The  desired  results  of  the  computation  are  the  mean  response  vec¬ 
tor,  |u{t)},  and  the  Covariance  Matrix,  Cov ({X(t)  j, {X{t) )  as  func¬ 
tions  of  time. 

No  attempt  is  made  to  model  the  distribution  of  the  response  vector 
or  to  detemine  the  correlation  between  values  at  different  times, 

3.2  TRANSFORMATION  OF  THE  EQUATIONS  OF  MOTION 

The  equation  of  motion  for  a  mass  excited  MDOF  system  is, 

[M]|X(t)|  +  ,C]{X(t)}  +  [Kj{X(t)  1  =  {F(t)  } 


(3-1) 
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where  [M]  =  System  Mass  Matrix 

[C]  =  System  Damping  Matrix 

[K]  =  System  Stiffness  Matrix 

{X(t)}  =  System  Displacement  Vector 
{F(t)}  =  Forcing  Function  Vector 

When  [l(]  and  [C]  are  deterministic  the  response  of  the  system  may 
be  obtained  by  several  methods.  For  a  deterministic  loading  {F(t)}, 
Equation  (3-1)  may  be  solved  by  a  variety  of  numerical  schemes  including 
modal  decomposition,  explicit  integration,  and  implicit  integration. 

When  {F(t)}  is  a  stochastic  process  modal  decomposition  is  the 
method  most  frequently  used.  However  when  [K]  and  [C]  are  matrices  of 
random  variables  the  solution  of  Equation  (3-1)  is  not  readily  avail¬ 
able. 

The  transformation  method  developed  in  Chapter  2  provides  a  new 
method  for  addressing  this  problem. 

Following  the  development  of  Section  2.2  results  in, 

{X(t)}  =  {u{c)}  +  {6(t)}  (3-2) 

Taking  the  expected  value  of  Equation  (3-2)  yields, 

E{X(t)}  =  E(y(t)}  +  E{6(t)}  (3-3a) 

=  {ii(t)}  (3-3b) 

To  develop  an  expresion  for  the  mean  square  value  of  the  response 
one  multiplies  Equation  (3-2)  by  its  transpose  yielding, 

{X(t)}{X(t)}^  =  {y(t)HM(t))'^  +  (M(t)}{6(t)}'^ 
t  {6(t)Hy(t)|'^  +  {6(t)}(6(t))^ 


(3-4) 
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Taking  the  expected  value  of  Equation  (3-4)  yields, 

E[(X(t)HX(t))T]  =  E[li.(t))|ii(t)}''']  +  E[{5{t)}(S(t)r]  (3-5) 

The  covariance  matrix  of  the  response  is  given  by, 

Cov[{X(t)},{X{t)}‘^]  =  E[{X(t)}{X(t)}'^]  -  E{X(t)}  .  E{X(t)}^  (3-6a) 

Cov[{X(t)},{X(t)}^]  =  E[{6(t)}{5(t)}''']  (3-6b) 

One  notes  that  the  expressions  in  Equation  (3-6b)  are  nxn  matrices 
for  a  system  with  n  degrees  of  freedom.  The  diagonal  terms  in  the 
matrices  are  the  variances  of  the  displacement  responses  at  the  degrees 
of  freedom,  and  the  off  diagonal  terms  are  the  covariances  between 
displacement  responses  at  the  degrees  of  freedom. 

A  similar  substitution  results  in  the  expression  for  the  mean  and 


covariance  terms  for  the  velocity  vector,  {X(t)},  which  are, 

{X(t)}  =  (u(t)}  +  {6(t)}  (3-7) 

E{X(t)}  =  (u(t)}  (3-8) 

Cov({X(t)},{X(t)}^)  =  E[{6(t)}{6(t)}^]  (3-9) 

The  expressions  for  the  mean  and  covariance  of  the  acceleration 
vector,  {X(t)}  are  similarly  derived  and  are, 

{X(t)}  =  (u(t)}  +  {6(t)}  (3-10) 

E{X(t)}  =  iu(t)}  (3-11) 

Cov({X(t)),{X(t)}^)  =  E({6(t)}(6(t)}'^)  (3.  12) 
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Representing  the  stiffness  matrix,  [K],  as  the  sum  of  [Kq]  and  [Kj^] 
with  [Kp]  deterministic  and  [Kp]  a  matrix  of  random  variables  with  mean 
zero  yields, 


[I'l  “  (i'n]  +  [l^ol 


(3-13) 


Taking  the  expected  value  of  Equation  (3-13)  yields 


E[K]  =  E[Kj3]  +  E[Kp] 


(3-14) 


And  since  [Kr]  is  mean  zero. 


E[K]  =  [Kj 


(3-15) 


To  develop  covariance  terms  one  must  consider  individual  elements 
in  the  matrix  [K].  Let  (K)^j  denote  the  element  in  the  i—  row  and  j — 
column  of  [ K] .  Then, 


(K),j  =  (K„),j  t  (K,),j 


(3-16) 


The  product  between  this  and  the  arbitrary  element  (K)j^ni  is. 


(3-17) 


Taking  the  expected  value  of  Equation  (3-17)  yields, 


E[(K)...(K)„  J  =  (K^).  (Kr,)„  +  E[(KJ.,(kJ,  ] 

'•  ij  '  Dmj  D  Jim  ^  R  ij  R  Jlm^ 


(3-18) 


because  the  random  components  are  mean  zero.  This  quantity  is  the 


correlation  between  the  random  variables  (K)..  and  (K)„  .  The  covari- 

1  j  im 

ance  between  random  variable  pairs  is  its  correlation  minus  the  product 
of  means.  In  the  present  case  this  is, 

Cov[(K),j,(K)^J  =  E[(K^),j(K,) J  (3-19) 

Similar  expressions  are  developed  for  the  damping  matrix,  [C],  with 
[Cq]  and  [C[^  ]  representing  the  deterministic  and  random,  mean  zero  terms 
respectively , 


[c]  =  [C;)]  +  [c,,]  (3-20) 

The  resulting  expression  for  the  mean  damping  matrix  is, 

E[C]  =  [C(,l  (3-21) 

The  covariance  between  terms  of  [C]  is, 

Cov[!C),j(C),J  =  E[(C„),j(C,)^l  (3-22) 

The  substitution  for  the  vector  of  forcing  functions,  {F(t)},  con¬ 
sists  of  the  sum  of  deterministic  mean  portion  {P(t)}  and  the  random 
portion  {q(t)].  Making  this  substitution  and  following  the  development 
above  results  in  the  following  expressions  for  |F(t)]  its  mean  and  co- 
vari ance. 

tF(t)}  =  {P{T)}  +  lq{t)}  (3-23) 

E{F(t)}  -  {P(T)}  (3-24) 

Cov(|F(tjHF(t)|^)  -  EI|M(t)K4(t)|'') 


(3-2b) 


71 


Substituting  the  results  of  Equations  (3-2),  (3-7),  (3-10),  (3-13), 
(3-20),  and  (3-23)  into  Equation  (3-1)  yields, 

=  {P(t)}  +  {q(t)}  (3-26) 

When  the  system  parameters  and  response  are  considered  correlated 
the  expected  value  of  Equation  (3-26)  is, 

[M]{y}  +  [Cpliil  +  [K(,]{y}  +  E([K^]{6})  +  E([C^]{6}) 

=  {P(t)}  (3-27) 

A  solution  for  {y(t)},  the  mean  response  vector  may  be  sought  oased  on 
this  expression. 

To  develop  an  expression  fvir  {6(t)}  one  may  subtract  Equation 
(3-27)  from  Equation  (3-26)  resulting  in, 

[Kllsl  ♦  [C„]li)  *  (Kj,](6}  1  ICRTl  *  {KRT)  = 

=  IrCO)  -  (KrIM  -  [CpJlil  (3-28a) 
where  {CRT|  •  [C^jl6l  -  E([C|,]{«))  (3-28b) 


{KRT|  =  [K^]{6l  -  E{[K^)l6l) 


(3-28c) 
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The  terms  in  Equations  (3-28b)  and  (3-28c)  may  be  considered  very 
small  and  may  be  neglected  under  the  arguments  presented  in  Sections 
2.2,  and  2.5,  resulting  in, 

[Mils}  +  [g{j)  t  (gjj)  =  {q(t)l  -  [K^Hul  -  [C|,l{;}  (3-29) 

The  forward  difference  method  is  used  to  solve  Equation  (3-27)  for 
the  mean  response  of  the  system.  The  expression  for  the  variance  of  the 
response  is  developed  by  applying  the  forward  difference  equations  to 
Equation  (3-*29),  yielding  an  expression  for  {6(t)}  in  terms  of  the 
underlying  random  variables,  and  {6(t)}  at  earlier  times.  This  expres¬ 
sion  is  then  mulliplied  by  its  transpose,  and  the  expected  value  of  the 
result  is  equal  to  th '  variance  at  time  t. 

The  details  of  this  development  are  provided  in  Section  3.4  after 
some  comments  on  the  development  of  the  system  matrices  in  Section  3.3. 
3.3  FORMULATION  OF  THE  SYSTEM  MATRICES 

The  finite  element  method  was  us<*d  to  develop  the  system's  stiff¬ 
ness,  damping  and  mass  matrices.  The  structure  of  interest  is  repre¬ 
sented  as  the  combination  of  a  number  of  place  frame  elements.  Each 
element  has  the  six  degrees  of  freedom  shown  in  Figure  3-1. 
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The  stiffness  tnatrix,  [K],  for  such  an  element  is. 


AE/L 

0 

0 

-AE/L 

0 

0 

0 

12EI/L^ 

6EI/l2 

0 

-12EI/L^ 

6EI/l2 

0 

6EI/l2 

4EI/L 

0 

-6EI/l2 

2EI/L 

a:/l 

0 

0 

AE/L 

0 

0 

0 

-12EI/L^ 

-6EI/l2 

0 

12EI/l2 

-6EI/l2 

0 

6EI/i2 

2EI/L 

0 

-6EI/l2 

4EI/L 

where  AE  =  Extensional  Rigidity  of  the  Element 
El  =  Flextural  Rigidity  of  the  Element 
L  =  Length  of  the  Element 

When  El  and  AE  are  random  variables  each  term  of  [k]  consists  of  a 
random  component  with  mean  zero  and  a  deterministic  component  equal  to 
the  values  obtained  from  Equation  (3-30)  evaluated  at  the  mean  values  of 
E I  and  AE.  Thus, 


[k]  =  [kRl  Mkol  (3-31) 

Takinc  the  expected  value  of  Equation  (3-30)  yields, 

E[kl  .  E(k„l  *  E(k„]  (3-32) 

And  since  [kj^]  is  '"ean  zero  and  [kj^]  is  deterministic, 


F  f  k  1  =  f  k 


(3-35) 
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When  individual  realizations  of  El  and  AE  appearing  in  Equation 
(3-30)  are  assumed  uncorrelated  then  the  elements  of  [kp]  are  indepen¬ 
dent  of  one  another.  This  assumption  is  made  for  ease  of  computation 
only.  In  fact  the  correlation  functions  between  the  terms  of  the  system 
stiffness  matrix  will  seldom  be  known  in  practice.  The  method  developed 
in  the  following  sections  may  be  applied  to  cases  where  the  terms  in  the 
stiffr-^ss  matrix  are  correlated,  provided  the  correlation  functions  are 
determined  and  stored  for  use  at  each  time  step.  This  becomes  apparent 
in  the  formulation  of  the  finite  difference  expressions  and  is  discussed 
in  detail  in  Section  3.4. 

The  system  stiffness  matrix,  [K],  is  formed  by  summing  the  contri¬ 
butions  of  each  element  for  each  degree  of  freedom.  When  the  orienta¬ 
tion  of  the  system  degrees  of  freedom  differ  by  an  angle  y  from  the 
element  degrees  of  freedom  the  matrix  defined  by  Equation  (3-30)  must  be 
rotated  to  the  proper  coordinates.  This  rotation  is  accomplished  by 
pre-  and  post-multiplying  [k]  by  a  transformation  matrix  [T].  When  one 
considers  an  element  oriented  as  in  Figure  3-2  the  stiffness  matrix 
expresses  the  relations  between  forces  and  displacements  in  the  X  and  Y 
di rection . 

The  transformation  matrix  [T],  used  to  develop  the  rotated  stif¬ 
fness  matrix,  [k'],  which  expresses  the  relations  between  forces  and 
displacements  in  the  X'  and  Y'  dif'ertions  is, 


[T] 


cosy 

-siny 

0 

0 


siny 

cosy 

0 

"t 

0 

1) 


0  0  0 

0  0  0 

1  0  0 

0  cosf  'iny 

9  -sin^  COSY 

0  0  0 


0 

0 

0 

0 

0 

1 


(3-34); 


,0 
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And  the  rotated  stiffness  r.iatrix  is. 


[k']  =  [Tf[k]CT]  (3-36) 

Since  [k]  has  a  deterministic  and  random  components  [k']  will  as  well. 
Thus  [k'l  may  be  expressed  as, 

[k'l  '  (Kr'!  +  [ko']  (3-36) 

with  the  random  component  [kp']  given  by, 

[kg']  =  [Tf  [KrKT]  (3-37) 

and  the  deterministic  component  [Kq*  ]  given  by 

[ko'l  MTl^[k[,](T]  (3-38) 

The  expected  value  of  [K']  is 

E[k'l  “  [ko'l  -  [Tl'|kQl[T]  (3-39) 

The  covariance  of  [k'][k']^  is  given  by, 

Cov[k'][k'f  =  t[k^'][k^'f  (3-40) 

Replacing  [kR*]  with  the  Equation  (3-27)  yields 

E[kR'|[kR')^  =  E(3|^kJ[Tl[Tf[k„|T[T] 


(3-41) 


/o 


However,  [T]  is  orthogonal  to  [T]^  so  [T][T]^  =  [I]  resulting  in, 

E[k«'][l<R'f  =  E[Tf[kR][ICRf[T]  (3-42) 


And  since  [t]  is  deterministic 

E[l<R'][kR’f  =  lTfE((k|,][k„f)[T]  (3-43) 


Resulting  in 


Cov[K'][k’f  =  [T]^{Cov[k][kf  )CT]  (3-44) 


However  due  to  the  manner  in  which  [Kj^  ]  appears  in  the  expressions 
developed  in  the  following  section  the  properties  of  interest  are  the 
variance  of  each  term  of  [K]  and  hence  each  term  of  [k'].  The  ij—  com¬ 
ponent  of  [k']  may  be  written  using  indicial  notation  as 


T.'^  j 

u  '^tp 


(3-45) 


Hence  the  Var(k^j)  may  be  expressed  as 

Var(k^)(Tpj)2  (3-46) 


when  the  terms  of  k  are  independent  as  assunied  earlier.  When  this  is 

not  the  case  the  values  of  the  covariances  of  fk.  k.  1  must  be  known. 

'■  1  j  tp 

When  iiie  correlation  between  A£  and  El  is  known  these  may  be  readily 
computed  from  Equation  (3-30).  However,  there  will  be  6**  of  them.  The 
final  asse-'bly  of  the  global  matrices  would  result  in  n^*  terms  of  the 
form  Cov(,Kijl<)^p}, 
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The  system  stiffness  matrix  is  assembled  by  summing  the  components 
of  the  element  stiffness  matrices  at  each  degree  of  freedom.  Thus  [K] 
may  be  written  as, 

CK]  =  I  [k'L  (3-47) 

i=l  ^ 

where  [k']^  is  the  stiffness  matrix  for  the  i— element,  and  z  is  the 
total  number  of  elements. 

Similar  expressions  for  the  random,  [Kj^  ]  and  deterministic,  [K^], 
components  of  [K]  are. 


And 

(3-49) 

i“l 

The  variance  of  each  term  K. .  is  the  sum  of  the  variances  of  k'  for 

ij  mp 

each  element  where  element  degrees  of  freedoms  m  and  p  correspond  to 
system  degrees  of  freedom  i  and  j  respectively. 

The  mass  matrix  for  a  frame  element  may  take  one  of  several  forms. 
The  one  chosen  for  use  here  is  presented  by  Cook  (33).  The  matrix  is 
diagonal  with  no  zero  entries  which  makes  it  trivial  to  invert  as  well 
as  offering  several  other  advantages  in  forming  the  finite  difference 
expressio.s.  The  terms  of  the  tnass  matrix  m  are  given  by, 

[mj  =  I'aiU  Vij  21U  210 
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where  M  is  the  total  mass  of  the  element.  When  required  (m)  may  be 
rotated  as  was  [k]  in  Equation  (3-35).  This  results  in, 

(m' )  =  [tT  (">)  CT]  (3-61) 

The  system  mass  matrix,  (M)  will  also  be  diagonal  and  is  formed  by  sum¬ 
ming  the  element  matrices  at  appropriate  degrees  of  freedom.  Thus  (M) 
may  be  expressed  as 

(M)  =  I  (m').  (3-52) 

i=l  ^ 


Where  (m)^  is  the  mass  matrix  of  the  i  element  and  i  is  the  total 
number  of  elements. 

The  system  damping  matrix  is  a  linear  combination  of  [K]  and  [M]. 
This  formulation  is  called  Rayleigh  damping  and  results  in, 

[C]  =  a[K]  +  0(M)  (3-53) 

where  a  and  6  are  deterministic  constants.  Substituting  Equations 
(3-20)  and  (3-13;  for  [C]  and  [K]  above  yields. 

[Clp  ♦  [Cj  =  a[K„l  +  a[KKl  *  6(M)  (3-54) 

Taking  the  expected  value  of  Equation  (3-54)  yields 

and  hence  the  random  component  of  [C]  is 


( 3-bb ) 
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Alpha  and  Beta  may  be  chosen  such  that  the  damping  of  the  system  is 
a  predetermined  portion  of  crilical  at  two  separate  frequencies. 
Considering  the  frequencies  uii  and  u)2  with  corresponding  percentages  of 
critical  damping  and  ei,  6  and  o  are  given  by, 


2  ( (1)2  £2  “  ^  ^ 

2  2 

012  “ 


(3-57) 


6  =  2ujjej^  -  (3-58) 

Once  the  system  matrices,  [Kq],  [Kj^],  [C^],  [Cj^],  and  (M)  have  been 
obtained  the  forward  difference  method  is  used  to  compute  the  mean 
(p(t)},  and  covariance  {X(t)  }{X(t)  of  the  response.  The  details  of 
this  computation  are  explained  in  the  following  section. 

3.4  APPLICATION  OF  FINITE  DIFFERENCES 

The  forward  difference  technique  is  used  to  express  the  equations 
of  motion  generated  in  Section  3.2.  The  resulting  expressions  are  thus 
evaluated  in  a  step-by-step  manner  to  yield  values  for  the  mean  and 
variance  of  the  response  over  times 

Let  tn  =  nAt,  r  =  0,  1,  ...,  At  a  "small  value",  be  a  discretized 
time  variable.  The  vector  forms  of  the  difference  expressions  developed 
in  Section  2,2  are, 

=  [f^l  (3-69) 

“-It  -  l*„))  (3-60) 

-  -J-  '  tf„!) 

At  ‘ 


(3-bl) 
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To  develop  the  relations  required  to  evaluate  {li(t^)}  and  its  deri¬ 
vatives  in  Equation  (3-27).  This  yields, 

(M)  ({v2t-2lvilMp„l)MCD]  (Ivil-Kl) 

+  =  {’’nl  '  -  £([<:„]{«„))  (3-62a) 

where  _  1 

(M)  =  ^^2  (M)  (3-62b) 

~  (3-62c) 

If  the  relationship  between  [Kj^j  and  [Cj^  ]  as  defined  in  Equation 
(3-56)  is  low  imposed,  Equation  {3-62a)  reduces  to, 

^^n+2'  ^  ^n^^  ^  ‘  ‘"n^^  ^  " 

=  {Pnl  -  ^[[KRlInn})  (3-63a) 

where  (n^^)  =  [6^}  +  o|(5^}  (3-63b) 

Solving  Equation  (3-63a)  for  yields, 

1p„*i1  =  (M)  ^  ((P„l  -  E([K„l|n„l)  ♦  (Aillw„^ll  [A^KrJ)  (3-64a) 

where  [A, ]  =  2  (H)  -  [0^]  (3-64b) 

[A,J  =  [C„|  -  (M)-[K„]  (3-<,4c) 


Replacing  [Cg]  with  Equation  (3-55)  yields. 


[Ajl  =  (2-StB)  (M)  -  (-ff)  [Kq]  (3-65a) 


[Aj]  =  (sAt-1)  (M)  +  (-fj.)  [Kp]  (3-65b) 


Equation  (3-64a)  provides  an  expression  by  which 
evaluated  in  a  step-by-step  manner  similar  to  that  used  for  the  SDOF 
system  discussed  in  Chapter  2.  The  values  for 
available  from  previous  time  steps. 

The  expected  value  of  may  be  developed  from  the  expres¬ 

sions  ^cr  {5p)  which  now  follows. 

The  process  {5(t)}  is  implicitly  defined  by  Equation  (3-29). 
Applying  the  forward  difference  relations  to  {6(t)}  and  its  derivatives 
results  in, 

m  (IVat  -  *  1«„!)  *  I^dKIViI  ■  t‘nl)  * 


-  [KrUzJ 


(3-66a) 


where 


|Z„1  =  ti-nl  ♦  “li-nl 


(3-66b) 


Solving  Equation  (3-65a)  for  (Sn+zl  yields, 


~  ^ ^'^n+1  (3-67) 


To  evaluate  the  expected  value  of  the  product 
writes  the  expression  for  t)y  evaluating  Equation  (3-67)  for  its 

kth  element  and  multiplying  by  ^nd  taking  the  expected  value. 

The  equation  for  tie  k^^  element  of  {6^+2}  is, 

^^n+2^k  "  ^^kk  ^^^n^k  "  ^^l^kr^VlV 


(3-68) 


Premultiolying  both  sides  by  [KR)ij  yields, 


^'^R^ij^^n+2^k  ""  ^^kk  j ^‘^n^k  "  ^“^R^ 


Taking  the  expected  value  of  Equation  (3-68)  yields, 


*  (\)k,E((^„),jU„,,),(  »  (*2)kpE((KR),j(«„)p))  (3-70; 


And  since  [Kr  ]  is  independent  of  (qn)  Equation  (3-70)  reduces  to, 


E((KR),,(6„,2)k)  - 


(''2>kp'«^Rln(«n)p') 


(3-71) 


Thus,  the  terms  of  the  E((K|^)^ evaluated  in  a 
step-by-step  manner.  However,  the  evaluation  will  require  n^  calcula¬ 
tions  for  each  of  the  n^  terms.  To  reduce  the  number  of  computations  a 
simplifying  assumption  may  be  imposed.  The  assumption  is  that  the 
deflection  at  degree-of -freedom,  (DOF),  k  is  independent  of  the  entries 
of  the  random  stiffness  matrix  at  dof's  other  than  k. 


This  may  be  written 


when  k  *  i 


or  k  j 


(3-72) 


When  k  =  i  the  E((K,^)^j(6^^.2^k^ 


*  (A2)kpE(«R)kj(«n)p)l 


(3-73) 


(no  sum  on  k  or  j) 


Due  to  the  independence  of  the  elements  of  [Kr ]  the  expected  value  of 
(KR)kj{KR)k£  is  zero  unless  it  =  j  then 


E((kp)kj(K„),,)  =  E(K„)7j 


(3-74) 


Furthermore,  by  taking  advantage  of  the  assumptu)n  expressed  in  Equation 
(3-71)  for  the  terms  (KR)ki (6n+l)r  and  KR)kj(6n)p, 

Equation  (3-73)  may  be  rewritien  as, 

^it^ikjlW’k)  =  f^)kk  [-Et^R)kj(E„)j  ^  (Al)kkE((SVj(Vl)k) 


'  ’'A3!k/i'.^'kji«„ij)i 


(3-7M 


(no  suii)  on  k  or  j)’ 


86 


When  we  consider  the  case  of  lc=j  Equation  (3-71)  may  be  written. 

Once  again  taking  advantage  of  the  independence  of  the  elements  of 
[K|^]  and  the  assumption  expressed  in  Equation  (3-72),  Equation  (3-76)  is 
reduced  to 

^  (Ai)jiE((KR),j(5„,i),)MA2)j^E((KR),.(6„)i)  (3-77) 

no  sum  on  j  or  i 

However,  the  expressions  in  Equations  (3-75)  and  (3-77)  allow  the 

Zn^  terms  of  (Ko)..(6  ).  for  k=j  and  k*i  to  be  evaluated  at  each  point 
K  1 J  n  K 

in  time.  When  i*j,  Equation  (3-71)  with  the  assumption  of  Equation 
(3-71)  yields. 


MA2),iE((K,),i(6n,i)il  (3-78) 

However,  if  one  considers  a  slightly  more  restrictive  assumption 
tnan  that  of  Equation  (3-72)  the  calculations  may  be  reduced  still 
further. 

Consider  the  product  coiiponent  of  the  ran¬ 

dom  portion  of  the  static  restoring  force  at  DOF  i  due  to  a  deflection 


87 


at  DOFj.  If  one  considers  the  product  (1^)1  clear 

physical  meaning.  One  may  therefore  say  that  the  terms  involving 

will  be  representative  of  the  uncertainties  in  actual  forces 
and  are  the  terms  of  primary  interest.  Since  the  terms  involving 
(K|^)ij(6n)i  have  no  such  clear  physical  meaning  it  is  assumed  they  will 
not  have  a  significant  effect  on  the  solution  of  the  problem.  The 
assumption  is  made  for  computational  convenience.  The  more  liberal 
assumption  of  j('^n)-i  1  ^nd  ^  equal  zero 

would  be  more  accurate  provided  computation  and  storage  capabilities 
exist  to  execute  it.  The  means  for  doing  so  are  provided  by  Equations 
(3-75)  and  (3-77).  The  full  set  of  n^  terms, 

be  assembled  at  each  time  step  by  Equation  (3-71),  provided  once  cgain 
tliat  sufficient  storage  and  computation  capabilities  are  available. 

Evaluating  Equation  (3-71)  under  the  constraint, 

E((K,)ij(5„)p  -  0  k.j  (3-79) 


yields. 


*  (Al)irEt(^),j(Vl)rl  *  jp^  >1  j(  «n)p  I  (3-30) 
Applying  Equation  (3-79)  to  the  Kr6  terms  on  the  right  hand  side 
of  Equation  (3-80)  yields. 


E((KR)ij(V2)j)  =  H((Kr,j)(Krj,))(ZJ^ 


t  (A,  )  .  .£( (K.),  (6  ^,)  )  +  (A.J  £((K^)  {6  )  1]  (3-Hl) 

'  rjj  ''  RMj'  n+l'j-'  ^  R'lj'  nM'*  '  ^ 


and  due  to  the  independence  of  the  terms  of  (Kr)  when  j^ti  this  yields, 


*  (Ai)jjE((K,),j(«„)j)]  (3-82) 

And  so  when  j('Sn)j)Q=  0*  the  case  of  zero  start.  Equation 

(3-82)  results  in. 


then  i=j  Equation  (3-81)  becomes. 


(3-83) 


(3-84) 


Hence  through  the  assumption  of  independent  elements  in  the  stiff¬ 
ness  matrix  and  the  constraint  of  Equation  (3-79)  the  expression  of 
Equation  (3-71)  involving  n^  terms  requiring  a  total  of  n^  (n^  per 
term),  computations  per  time  step  is  reduced  to  n  terms  involving  a 
total  of  n  computat’ons  per  time  step. 

Equation  (3-34)  now  allows  one  to  evaluate  Equation  (3-64a)  for  the 
mean  response  since  is  a  vector  whose  i— term  is  given  by 


‘^**'R*H**n+l 


it) 


(E[RRl|n„}),  = 


At 


(3-85) 
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To  develop  an  expression  for  the  variance  of  the  response  one  first 
multiplies  Equation  (3-67),  the  expression  for  transpose 

and  then  takes  the  expected  value  of  the  resulting  expression.  This 
results  in, 

-  [Ejf  -  [EJ  t  [E^] 

+  [Ejf  *  [E3I*  [Ejf  *E[K^]tZ„Hz/l'-Rf 

-  [E4]  -  [E  ,f  -  [E5I  -  [Ejf  *  .  [E.ft  [Eg] 


♦  E[A2)[S„H«„r[A2])  (M)-'  (3-864) 

where  [EJ  -  E(|q„}{Zj^[K|,]^  (3-86b) 

[Ej]  '  E(lq„}{sJ^[A2f )  (3-86C) 

[E3I  -  Eilq„l|5„^ll^lAif)  (3-86d) 

[E4I  '  E([K„HZ„){S„^j}^[Aj]^)  (3-86e) 

IE5!  •  E([h^]|Zj|S„l^[A2]^)  (3-86f) 

[Eel  •  E([A,l{6„^;)|6„lT[A2f) 


(3-86g) 


9 


■M 

Since  the  loading  and  stiffness  are  independent  [Ej]  may  be 

written, 

V- 

(Ell  =  Elq„)lZ„lV,K„]'^ 

(3-87) 

And  since  {q^}  and  [K^^]  are  mean  zero. 

[Ej]  =  [0] 

(3-88) 

Upon  examination  of  Equation  (3-67)  one  notes  that  {6}p  is  independent 

of  {q^  and  Hence  [E^]  and  [E^]  may  be  written  as 

[Ejl  =  E|q„)El6„lT[A2]’^ 

(3-89) 

9 

[E3]  •  E(q„)E{6„,;lV,l^ 

(3-90) 

^  ■' 

And  since  {q^}.  {"5^}  und  are  all  mean  zero, 

[E2]  =  [0] 

(3-91) 

j» 

and 

[E3]  =  [0] 

(3-92) 

• 

To  evaluate  [E^]  one  may  write  out  the  expression  for  its  in 

term, 

ijth 

(Eqlij  ■  E((KR)ii(Z„),(Vl>p(V)pj) 

(3-93) 

0 

• 

«• 

» 

i'"  J'" i 


'r*  v  ■  T."  *  u:^ U^P  ii  i  ■  1',^ 
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However,  from  Equation  (3-84)  this  reduces  co, 

(E4)iJ  =  E((K,)„(V))Z„/iiJ  (3-94) 

Simi  larly, 

=E«KR)„(Vl)Z„i''2'j  (3-36) 

The  matrix  [Eg]  may  be  written  as 

[Efi]  =  [A,]EH«„^lH«„l^)[A2f  (3-96) 

Note  that  {fin+ll  ’’’^y  be  written 

IVil  “  *  [Alltel  *  [AslIVil^  (3-33) 

Postmultiplying  the  above  by  {6^}^  yields, 

lViH4„r  ^  -  t^Rl|Z„.lH4„}'  *  IA2H4„H«„}' 

*  [A2H«„.iII«„1^)  (3-98) 

Taking  the  expected  value  yields, 

E|«„*iH«„r  =  (Mr‘[0  -  ECKp]|Z„.jH4„l^  *  [Al)Ei6„H«„l' 

♦A2(e|«„HVi>')'] 


(3-99) 
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The  expression  is  merely  the  transpose  of  the  left 

T 

nand  side  evaluated  at  the  previous  time  step.  The  E  [ [Kj^  ]  {6^} ’ } 
is  a  diagonal  matrix  whose  entries  are. 


=  (Z„.l)iE((KR)(«„)i)  (3-100) 

The  expression  E [Kp may  be  written  in  tensor  nota¬ 


tion  as 


E((KR)i,{Z„),(Z„)p(KT)pjl  =  (Z„),(Z„)pE(K„,4j)  (3-101) 


J 

And  since  E(K[^^  =  0  unless  j-i , 

diagonal  matrix  with  entries  given  by 


p=t,  Equation  (3-101)  reduces  to  a 


T  _ 


2 

E(Ko.  .)  Z 
'  Rij'  nj 


Equation  (3-86c)  may  now  be  evaluated  at  each  time  step  yielding  an  • 
expression  for  the  Cov((X(t)  }{X(t)  p)  at  euCh  time  step. 

The  numerical  algorithm  for  performing  the  computation  is  presented 
in  Section  3.5.  The  computation  provides  the  expected  value  of  the  dis¬ 
placement  response  at  each  degree  of  freedom  and  the  covariance  between 
the  responses  at  various  degrees  of  freedom. 
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3.5  THE  MULTIPLE  DEGREE  OF  FREEDOM  COMPUTATIONAL  ALGORITHM 

The  following  algorithm  is  used  to  evaluate  the  mean  response, 
iy(t)},  and  the  covariance  of  {X(t)}  using  the  forward  difference  equa¬ 
tions  developed  in  Section  3.4,  The  algorithm  is  implemented  in  the 
code  RMDOF.  A  listing  of  the  code  is  provided  in  Appendix  B. 

The  algorithm  is  based  on  a  zero  start  condition  ind  proceeds  as 
fol lows. 


1.  Read  the  problem  description  parameters,  the  number' 

of  nodes,  the  number  of  elements  and  the  number  of  restraints. 

2.  Read  the  nodal  locations. 


3.  For  each  element,  read  the  nodes  to  which  it  is 
connected,  the  material  properties,  and  their 
uncertainties.  Formulate  the  element  stiffness 
matrix  and  compute  the  variance  of  each  term  in  the 
stiffness  matrix.  Sum  these  into  the  global 
matrices.  Formulate  the  element  mass  inotrix  and 
sum  the  results  into  the  global  mass  matrix. 

4.  Read  the  damping  parameters,  -iij,  S2  calculate  a 

and  8  from  Equations  (3-57)  and  (j-BS)”. 

5.  Formulate  the  problem  solution  natrices  [A^j  and  [Aj?] 
defined  in  Equations  (3-65a)  and  (3-55h). 


6.  Impose  the  constraint  conditions. 

For  each  time  step  do  the  followiruj: 

a)  Update  time  and  the  index  on  the  arrays  {un+i]» 

(Pnl.  Ef  1 > 

l^n+xi  *!^n}^] »  E[  ) , 

El[KRl{5n^l}j. 

b)  Form  E( 6^+21)  >  rom  Lquatiofi  ('3-85'). 

c)  Compute  L[  ( Sf.-l } l-nr  1  from  FMuation  (i-93). 


7. 
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d)  Compute  {un+2}  Equation  (3-64a). 

e)  Compute  E [ {6^+2 }f 6^+2}^]  from  Equation  (3-86a). 

f)  Write  appropriate  results. 

8.  Repeat  7  until  the  problem  is  completed. 

3.6  NUMERICAL  EXAMPLES 

This  section  presents  the  results  of  two  example  problems  which 
illustrate  the  aril)sis  of  MDOF  systems  with  random  stiffness  and  damp¬ 
ing.  In  both  examples  the  random  portion  of  any  element  of  the  global 
stiffness  matrix,  [K],  was  assumed  to  be  perfectly  worrelated  with  the 
random  portion  of  the  corresponding  element  in  the  damping  matrix. 

This  assumption  of  independent  terms  has  a  serious  shortcoming  in 
that  a  rigid  body  displacement  of  any  element  will  contribute  to  the 
variances  of  the  displacements. 

The  random  properties  of  each  element  of  the  stiffness  and  damping 
matrices  were  assumed  independent  of  all  other  entries  of  these  mat¬ 
rices. 

These  assumptions  were  made  for  reasons  of  computational  covenienr.e 
only.  When  one  is  modeling  simple  structures  such  as  homogeneous  bars 
and  beams  where  the  random  properties  are  assumed  nighly  correlated 
along  their  length,  these  assumptions  may  not  be  valid. 

However,  in  the  examples  below  the  elements  are  considered  to  be 
condensed  representations  of  many  beams,  columns  and  girders.  For 
example,  ''f  one  were  to  model  a  large  building  the  axial  and  flextural 
properties  of  several  bays  and  stories  may  be  combined  and  represented 
by  a  single  frame  element.  Thus  a  large  structure  may  be  modeled  in  an 
economical  manner.  Sucn  subassembly  modeling  is  common  fur  lar-ye  struc¬ 
tures  including  aircraft,  snips,  and  dams  as  well  as  buildings. 


1 
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EXAMPLE  1  -  Deflection  of  a  Non-Prismatic  Cantilever  Beam, 

The  structure  considered  here  is  shown  in  Figure  3-3.  All  dimen¬ 
sions  are  unitless.  The  structure  possesses  six  degrees  of  freedom. 
However,  the  discussion  here  will  concern  the  midpoint  and  tip  deflec¬ 
tion.  The  mean  and  standard  deviation  of  the  tip  deflection  are  plotted 
in  Figure  3-4.  The  mean  and  standard  deviation  of  the  midpoint  deflec¬ 
tion  are  plotted  in  Figure  3-5. 

The  correlation  coefficient  between  the  two  is  plotted  in  Figure  3-6. 

The  mean  and  standard  deviations  for  this  problem  follow  the  same 
basic  trend  as  that  observed  for  the  SDOF  examples,  namely  that  the 
standard  deviation  curve  has  the  same  basic  shape  as  the  mean  curve. 
However,  the  standard  deviation  curves  do  not  lag  behind  the  mean  curves 
as  in  the  SDOF  cases.  Furthermore  the  standard  deviation  does  not  dis¬ 
play  the  secondary  peaks  which  were  indicated  in  the  SDOF  case. 

The  ratio  betweei  the  midpoint  and  tip  deflections  appear  to  be 
about  a  factor  of  1/8.  The  standard  deviation  at  the  tip  is  approxi¬ 
mately  1/3  the  value  of  the  mea/.  The  standard  deviation  at  the  mid¬ 
point  on  the  other  hand  is  approximately  1/2  the  mean.  The  coefficient 
of  variation  of  the  stiffness  and  damping  were  0,15.  Since  more  elem¬ 
ents  contributed  to  the  stiffness  at  the  midpoint,  and  the  stiffness  of 
the  elements  is  Independent,  it  is  not  surprising  that  the  midpoint  has 
a  larger  coefficient  of  variation. 

The  coefficient  of  variation  at  the  tip,  however,  is  still  approxi¬ 
mately  double,  (1/3  vs  0.15}  what  would  be  expected  based  on  the  results 
of  the  SDOF  computations.  However  the.  assumption  of  independence  of 
terms  within  the  stiffness  matrix,  means  only  positive  terms,  speci¬ 
fically  those  of  Equation  (3-102':  enter  the  computations.  If  the  terms 


of  the  stiffness  matrix  were  to  be  considered  perfectly  correlated  then 
negative  terms  would  enter  the  computation  as  well,  reducing  the  overall 
variance.  The  SDOF  approximation  to  this  same  problem  would  in  effect 
assume  perfect  correlation  between  the  random  variables,  thus  the 
results  would  not  be  comparable. 

The  correlation  between  the  tip  and  midpint  deflections  approaches 
a  constant  value  of  0.615.  The  correlation  remains  positive  since  the 
deflections  always  are  in  phase,  and  tied  to  one  another  in  that  the 
deflected  shape  of  the  structure  is  well  defined  for  this  particular 


EXAMPLE  2  -  Axial  Response  of  a  Prismatic  Bar. 

The  structure  of  interest  for  this  problem  is  shown  in  Figure  3-7. 
All  dimensions  are  unitless.  The  response  parameters  of  interest  were 
the  midpoint  and  tip  deflections.  Two  loading  conditions  were  investi¬ 
gated,  the  first  being  a  loading  applied  suddenly  to  the  midpoint  and 
the  second  being  a  loading  applied  suddenly  to  the  end.  The  loading  in 
both  cases  had  a  value  of  .1. 

For  the  case  of  the  loading  at  the  midpoint  the  mean  and  standard 
deviation  for  the  deflection  at  the  tip  are  plotted  in  Figure  3-8.  The 
mean  and  standard  deviation  for  the  midpoint  are  plotted  in  Figure  3-9, 
and  the  correlation  between  them  is  plotted  in  Figure  3-10.  The  same 
information  for  loading  at  the  tip  is  plotted  in  Figures  3-11,  3-12,  and 
3-13. 

Examination  of  Figures  3-8  and  3-9  indicate  that  the  standard  devi¬ 
ation  of  the  tip  deflection  is  significantly  higher  than  that  at  the 
midpoint  even  though  the  deflection  are  nearing  their  final  values.  The 
increase  in  the  tip  deflection's  standard  deviation  is  due  to  the 
assumption  of  independent  entries  of  the  stiffness  matrix.  Hence  while 
the  end  element  has  undergone  only  a  rigid  body  deformation  once  the 
oscillations  are  damped  out  this  has  still  contributed  to  the  final 
standard  deviations. 

The  resulting  means  are  oscillating  about  a  point  which  is  approxi¬ 
mately  7%  greater  than  their  deterministic  static  deflection.  This  is 
due  to  the  impact  of  the  uncertainties  on  the  mean  response. 

The  correlation  coefficient  for  both  cases  is  approaching  a  cons¬ 
tant  value  of  0.5.  The  standard  deviations  for  the  case  of  the  end 


gure  i-iJ  Correlation  Between  the  Midpoint  and  En 
Bar  With  the  Loading  Applied  at  the  End 
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loading  are  higher  than  those  for  the  midpoint  loading.  This  is  due  to 
larger  deflections  on  the  second  half  of  the  bar  for  +he  case  of  an  end 
loading.  The  mean  value  of  the  midpoint  deflection  is  also  greater  for 
the  case  of  an  end  loading  than  the  case  of  a  midpoint  loading,  but  by 
only  a  few  percent.  This  is  due  once  again  to  the  effect  of  the  ran¬ 
domness  of  the  system  parameters  on  the  mean  values. 


iiijiinM 
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CHAPTER  4 

PASSAGE  PROBABILITIES  FOR  NONSTATIONARY  RANDOM  PROCESSES 

4.1  INTRODUCTION 

This  chapter  presents  a  new  method  for  estimating  the  probability 
of  passage  beyond  a  barrier  level,  h,  by  the  nonstationary  random 
process  X(t). 

When  the  barrier  b  represents  a  predetermined  failure  level  then 
the  probability  of  passing  beyond  this  barrier  is  identical  to  the 
probability  of  failure. 

The  method  of  computing  the  passage  probabilities  for  the 
nonstationary  random  process,  X(t),  presented  in  this  chapter  is  a  new 
extension  of  existing  methods  used  for  stationary  processes.  The 
extension  involves  the  application  of  conditional  probabilities  and  the 
assumption  that  the  process  X(t)  appears  stationary  during  short  time 
intervals. 

The  calculation  of  passage  probabilities  for  stationary  processes 
is  discussed  in  Section  4,2.  The  extension  of  this  to  nonstationary 
processes  is  presented  in  Section  ^,3  with  the  computational  algorithm 
for  implementing  the  new  method  presented  in  Section  4.4.  Numerical 
examples  of  the  application  of  the  technique  are  provided  in  Section 
4,5. 

4.?  PASSAGE  PHOBABILiriES  F0»  STATIONARY  RANDOM  PROCESSES 

The,  methodology  for  approximating  the  passagt  probabilities  for 
stationarv  random  processes  has  been  well  developed  by  Rice  (21)  and 
othr-rc.  Mio  metnod  developed  by  Rice  is  summarized  below. 


Ill 


If  one  considers  the  random  process  X(t)  to  be  weakly  stationary 
then  its  mean,  vij^(t),  and  variance,  are  constant  values  for  all 

values  of  t.  Furthermore  the  process  defining  its  derivative  with 
regard  to  t,  X(t),  is  independent  of  X(t). 

The  characteristic  of  interest  of  X(t)  is  the  collection  of  its 
random  peaks  which  fall  above  a  predetermined  level,  b,  where  b>y  .  Let 

A 

N(b;T)  be  a  counting  process  associated  with  the  crossing  of  the  barrier 
level  b.  Thus  every  crossing  of  b  by  X{t)  is  counted  by  N(b,t).  Speci¬ 
fically  it  is  assumed  that  N(b,t)  is  a  Poisson  counting  process. 

The  events  counted  by  a  Poisson  counting  process,  in  this  case  the 
crossings  of  b  by  X(t),  are  assumed  to  have  the  following  properties: 

a)  The  events  occur  independently,  hence  knowing  an  event  has  just 
occurred  provides  no  information  regarding  the  probability  of 
when  the  next  event  will  occur, 

b)  The  probability  of  one  event  occurring  in  a  time  interval  dt  is 
X5(t)dt  where  dt  <<  1,  and  X(t)b  is  the  mean  crossing  rate. 

When  X(t)  is  stationary  X^Ct)  =  Xjj  a  constant. 

c)  The  probability  of  more  than  one  event  occurring  in  a  time 
interval  dt  is  negligible. 

Assumption  (a)  above  is  arbitrary  for  structural  response  since  the 
realizations  of  the  process  X(t)  are  not  independent.  When  X(t)  is  at 
or  near  a  barrier  level  then  X(t+dt)  will  be  as  well  for  dt  <<  1.  This 
effect  is  called  "clumping"  and  is  addressed  by  Vanmarke  (4)  and 
others.  However,  this  effect  was  neglected  by  Rice  and  will  be  neglec¬ 
ted  in  the  following  sections  as  well. 
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The  Poisson  process,  N(b,t),  has  the  following  probability 
function, 

(Xbt)^ 

Pj^(j,t)  =  exp  (-Xj^t)  —  (4-1) 

where  j  is  the  number  of  events  counted,  that  is  the  number  of  times 
X(t)  crosses  the  barrier  level  b. 

The  probability  of  not  crossing  the  barrier  in  the  interval  (o,t) 
is  provided  by  Equation  (4-1)  with  j=o.  This  yields, 

P|^(o,t)  =  exp  (-X^t)  (4-2) 

Let  P|^(F,t)  be  the  probability  of  crossing  the  barrier  level  at 
least  once  in  the  increment  (o,t),  then, 

P,^(F,t)  *  1  -  exp  (-X^t)  (4-3) 


All  that  now  remains  is  to  compute  the  mean  crossing  rate,  X^. 
The  mean  crossing  rate  at  any  time  (t),  Xb(t),  if  given  by. 


x^(t)  =  /  i 


p^-^  (b,x,t)dx 


(4-4) 


where  p„,,  (b,x,t)  U  the  joint  probability  density  function  of  X(t)  and 

A  A 

X(t)  evaluated  at  X(t)  =  b.  Since  the  random  processes  X(t)  and  X(t) 
have  been  assumed  stationary,  Xt5(t)  =  X^,  and  Equation  (4-4)  may  be 
reduced  to. 


X 


b 


p^^  (b.xldx. 


(4-5) 


I 


y  -  ^  V%  ■  *ri  •?'  ^  w^*  v'l  V"  !?F  T*tr'«'7«  L'"" 
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When  X(t),  X(t)  are  independent,  mean  zero,  normal  random  processes 
their  joint  density  is  given  by. 


''s  •  "  W^a~  ^~2  ■  "S’  ^ 

XX  XX  2a  2a 


(4-6) 


Substituting  this  expressioti  into  Equation  (4-5)  yields, 


X 


b 


(4-7) 


Equation  (4-7)  now  defines  the  mean  crossing  rate  of  a  barrier  b  by 
the  stationary  mean  zero  random  process  X(t).  However,  this  rate 
reflects  crossing  the  barrier  from  bcth  above  and  below  and  since  the 
engineer  is  only  concerned  with  crossings  from  below  the  crossing  rate 
must  be  modified.  Fortunately  this  is  straightforward  since  for  a 
stationary  random  process  all  crossings  from  below  must  be  paired  with  a 
crossing  from  above.  Hence  the  mean  c^nssing  rate  from  below,  x"^  is 
merely  half  the  total  mean  crossing  rate. 


X 


+ 

b 


(4-8) 


The  resulting  probability  of  crossing  the  barrier  level  b  from 
below  within  the  increment  (o,t)  is  then 


P||j  (F,t)  =  1  -  exp  (-xjt) 


(4-9) 


These  results  may  now  be  expanded  to  include  the  present  case  of 
interest,  the  no-^  ^tati onary  response  of  a  structure  with  random 
stittness  and  dampiriq. 


4.3  THE  NONSTATIONARY  CASE 


A  new  method  for  computing  the  probability  of  having  crossed  a 
barrier,  b,  at  least  once  in  the  interval  {o,t)  for  a  nonstationary 
process  will  now  be  presented.  The  barrier  b  is  greater  than  the 
largest  values  of  When  Uj^(t)  is  greater  than  b,  the  probability 

of  passage  is  assumed  to  be  1. 

The  method  is  formulated  under  the  following  assumptions, 

a)  The  nonstationary  random  process  X(t)  and  X(t)  are  continuous. 
Hence  their  values  may  not  instantaneously  jump  at  any  time. 

This  is  completely  appropriate  for  processes  defining  stru  ‘ural 
displacement  and  velocity  response  under  most  loadings  of 
interest,  but  excludes  processes  such  as  Brownian  motion  and 
quite  possibly  buckling  responses. 

b)  The  means  and  variances  of  X(t)  and  X(t)  vary  only  slightly  over 
increments  of  time.  At,  which  are  small  in  comparison  to  the 
response  time  of  the  structure. 

c)  X(t)  and  X(t)  are  correlated  normal  random  processes  whose 

means,  standard  deviations,  0j|(t).and 

correlation  coefficient  r{t)  are  known. 


Such  a  function  obviously  does  not  allow  the  direct  use  of  Equation 
(4-9)  to  calculate  the  probability  of  passage  beyond  a  barrier,  b. 
Equation  (4-4)  will  still  provide  the  mean  crossing  rate  at  any  time  t. 
However,  consider  tne  behavior  of  X(t)  over  a  short  period  of  time.  For 
example  the  interval  (T,T+At).  Then  as  At  becomes  small. 


1.^(t)  =  Uj^(T  +  At) 

(4-10) 

and 

(4-11) 

Considering  the  approximations  of  Equations  (4-10)  and  (4-11)  the 
mean  and  standard  deviation  of  X(t)  within  the  increment  (T,T+At)  may  be 
approximated  by, 


Pj^(t)  =  yj^(T)  T  <  t  <  T+At 

(4-12) 

+  i.  ^  £  T+At 

(4-13) 

The  values  for  pj^(T+At)  and  Oj^(T+At)  could  have  been  used  for  the 
right  hand  sides  of  Equations  (4-12)  and  (4-13).  However,  this  would 
only  result  in  the  final  results  being  shifted  to  the  left  on  the  time 
axis  by  an  amount  At. 

The  similar  expressions  for  the  process  X(t)  are, 

Ux(t)  “  Ux(t)  t  £  t  £  T+At  (4-14) 

and 

aj^(t)  =«  aj^(T)  T  £  t  <_  T+At  (4-15) 

The  approximations  of  Equations  (4-12)  through  (4-15)  will  have  the 
largest  error  when  the  response  begins,  that  is  as  the  structure  moves 
away  from  the  deterministic  zero  start  condition.  However,  since  the 
probability  of  exceeding  a  barrier  in  a  time  At  when  starting  from  a 
zero  position  is  approximately  zero,  the  errors  occurring  at  early  times 
will  in  effect  be  "laundered  out"  by  the  computations. 

Consider  the  random  processes  X(t)  and  X(t),  with  the  properties 
outlined  above.  The  mean  upcrossing  rate  of  a  barrier  b  at  an  instant 


where  ^(b,x)dx  is. 


b-r  2 


Mv 


XX 


X  i-  ■  2-'  (V^)(- 

A,T  ^"^^X^Xp  2p2  '^x  '^X 


x-pO 


^  (V)  n 

°x 


and 

PX  = 

Px  =  yx(T) 

ax  =  a^{^) 

r  =  r(T) 

Letting, 

b-Py 

b'  =  ^ 

pax 

and. 

2  _ 

Equation  (4-16)  may  be  written. 


xj(r)  /"  (pOj^Z+p-)  exp{-^  [z2.2rb'Z+b‘2]}dZ 


Completing  the  square  of  the  exponential  term  in  the  integral  yields. 


Xb(-)  = 


exp  {-  Y  (p^b'^1} 

2TTa  , 

X  -Wjj 


00 

/  (pa^z  +  exp  {-  Y  (Z-rb')^}dz  (4-21) 


'3 


Now  letting, 


Z-rb*  =  S 


(4-22) 


^  -  rb'  =  S, 


(4-23) 


where  $[_  is  the  new  lower  bound  of  integration  yields. 


^  exp{ -  4  (p^b'^]}  «  , 

Xb(T)  =  - ^ -  /  (po^(rb+s)+yj(  exp{- [S)^}dS  (4-24) 


Expanding  this  expression  to  two  integrals, 

exp{-4p^b'^}  «  , 

^b(")  = - -  [(pOxPb+w;)  I  exp{-  2  S  }dS 


+  pojl  f  S  exp{  -  Y  5^} 


(4-25) 


Upon  solution  this  yields. 


^  exp( -  4  P^b'^}  , 

- lit -  ["“x  «pf-  2  1  -  s^f 


-  S.pot  /  expj-  4  S^}dS]  (4-26) 


The  integral  on  the  right  hand  side  of  Equation  (4-24)  may  be  evaluated 
in  terms  of  the  error  function,  erf{x)  where, 


1  X  1 

erf(x)  =  . . .  /  exp(-7y^)dy 

^TT  0 


(4-27 


erf  (»)  =  j 


(4-28 


Hence,  Equation  (4-24)  becomes. 


pa.  exp  {-  i  p^b'^}  , 


(4-29a 


where. 


G=-^  -  erf(S|^)  o<S^ 


(4-29b 


G  =  ^  +  erf(SL) 


S|^  <  0 


(4.29c 


Now  replacing  b'  with  its  definition.  Equation  (4-18),  Equation  (4-29a) 
becomes 

^b<'l  =  ^  !-  7  1-  7  M  1^-30 

where  $[_  is  defined  by  Equation  (4-23)  and  G  is  defined  by  Equations 
(4-29b)  and  {^-29c).  The  value  of  erf(x)  is  available  as  a  computer 
generated  function. 


Then  the  probability  of  crossing  the  barrier  b  at  least  once  during 
the  interval  (x.T+At),  P(b,T,At),  is 

P(b,T,At)  =  1  -  xJ(T)At  .  (4-31) 

Having  the  probability  of  at  least  one  crossing  occurring  in  the 
time  (o,t)  may  now  be  developed.  Let  t^=nAt,  where  n  is  a  positive  in¬ 
teger,  and  At  is  a  small  positive  time.  Let  the  event  a  be  the  crossing 
of  the  barrier  at  lest  once  in  the  interval  (o,tn).  Let  the  event  3 
be  the  crossing  of  the  barrier  at  least  once  in  the  interval  (o.tp.i). 

The  probability  of  a,  P(a)  may  be  written  as  a  the  sum  of  probabil- 
ilities  conditioned  on  3.  This  results  in, 

P(a)  =  P(a/3)P(3)  +  P{ct/3‘^)P(3‘^)  (4-32) 

where  P(a/3)  is  the  probability  that  a  occurs  given  that  3  has  occurred 
and  P(a/3^)  is  the  probability  that  a  occurs  given  3  did  not  occur. 

Since  P(a/3)  is  the  probability  of  at  least  one  crossing  in  the 
interval  (o,tn)  given  at  least  one  crossing  in  (o,tn-i), 


And  Equation  (4-32)  becomes. 


P(«)  =  P(3)  +  P(cx/3‘^)P(e'^) 


(4-34 


However  since  the  only  way  for  a  passage  to  occur  in  the  interval 

(o,t^)  given  that  none  occurred  in  the  interval  (o,t^_^)  is  for  the  pas 

sage  to  occur  in  the  interval  (t^  given  by  Equation 

(4-31)  evaluated  for  T=t  Hence 

n-1 


P(a)  =  P(3)  +  (l-xJ(t^_j)At)P(3^) 


(4-35 


Letting  PYi,(j)  be  the  probability  of  X(t)  achieving  at  leas^  one 


crossing  of  the  barrier  b  by  the  end  of  the  j— -  time  step  yields. 


P(a)  =  Pxb(") 


P(3)  =  Pxb(n-l) 


P(3^)  =  1-P(3)  =  l-Pxb(^-l) 


(4-36a 


(4-36b 


(4-36C 


Substituting  the  above  into  Equation  (4-35)  yields  the  recurrence 


relation, 


Pxb(n)  =  ^  (^+S^t^_^)At)[l-P^^(n-l))  (4-37 


This  expression  may  now  be  evaluated  in  a  step-by-step  manner  to 

yield  the  probability  of  X(t)  crossing  a  barrier  b  by  the  time  tn. 

The  computational  algorithm  for  implementing  Equation  (4-37)  is 

provided  in  Section  4.4  and  some  numerical  examples  are  provided  and 
discussed  in  Section  4,4. 
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4.4  THE  COMPUTATIONAL  ALGORITHM 

The  following  algorithm  was  developed  to  evaluate  Equation  (4-37) 
to  compute  the  probability  of  passage  beyond  a  barrier  b.  The  algorithm 
is  implemented  in  the  code  RSDOF.  A  listing  of  RSDOF  is  included  in 
Appendix  A. 

1)  Let  n=o,  set  barrier  level  of  interest 

2)  Let  n=n+l,  tn=nAt 

3)  Get  the  values  for 

r(t^)  as  computed  in  Chapter  2. 

4)  Evaluate  B  of  Equation  (4-20c)  to  determine  wh’cn 
expression  for  to  use 

5)  If  B  <  0  use  x^(t^)  as  defined  in  Equation  {4-30) 

6)  If  B  >  0  use  X^(t^)  as  defined  in  Equation  (4-29) 

7)  Compute  the  probability  of  crossing  a  barrier  b  during 
this  time  step  by  Equatica  (4-31). 

B]  Compute  the  probability  of  crossing  the  barrier  in  the 
interval  (o,tp)  from  Equation  (4-37) 

9)  Go  to  Step  2  until  problem  complete. 

The  results  obtained  by  implementing  the  above  algof'ithm  are 
discussed  in  the  followng  sections. 


4.5  NUMERICAL  EXAMPLES 


This  section  presents  the  results  of  five  calculations  using  the 
method  developed  in  the  previous  sections  of  this  chapter. 

The  first  example  consisted  of  three  calculations.  The  probability 
of  having  achieved  at  least  one  crossing  of  a  barrier  was  computed  for 
Cases  7,  8,  and  9  of  Example  2  in  Chapter  2.  These  cases  had  the  param¬ 
eters  shown  in  Table  2-2.  The  barrier  level  was  0.36,  which  is  90%  of 
the  peak  response  for  the  undamped  system.  The  only  variation  between 
cases  was  the  correlation  between  the  random  stiffness  and  damping. 

The  resulting  passage  probabilities  are  plotted  in  Figure  4-1.  The 
prob'jbi lities  rapidly  reach  a  value  that  is  close  to  their  final  value 
as  the  response  reaches  its  first  peak.  The  probabilities  increase 
slightly  with  each  additional  peak  but  are  assymt'tical ly  approaching  a 
final  value. 

The  probability  of  passage  for  the  case  of  a  negative  correlation 
between  the  stiffness  and  damping  is  the  largest  of  the  three  cases  fol¬ 
lowed  by  the  case  of  independent  stiffness  and  damping,  with  positively 
correlated  stiffness  and  damping  yielding  the  smallest  probability  of 
passage. 

The  second  example  is  the  probability  of  the  system  define!  in 
Example  3  of  Chapter  2  exceeding  the  peak  level  obtained  by  the  deter¬ 
ministic  calculation. 

The  resulting  passage  probability  is  plotted  in  Figure  4-2.  These 
values  are  for  the  system  with  8%  damping  and  a  coefficient  of  variation 
for  the  stiffness  of  0.15,  the  values  used  in  Chapter  2. 
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Hence  there  is  almost  a  40%  chance  of  this  structure  failing  when 
failure  is  defined  as  exceeding  the  peak  determinist  response.  The 
second  calculation  computed  the  probability  of  exceeding  the  value  of 
0.2884",  which  is  the  deflection  resulting  from  twice  the  peak  load  when 
it  is  applied  as  a  static  load. 


^  2  X  300,000 

2.08  X  10^ 


(4-38' 


The  results  of  this  computation  are  plotted  in  Figure  4-3.  The 
shape  of  tne  curve  is  similar  to  that  of  Figure  4-2;  however,  the 
resulting  probability  of  passage  is  now  only  about  .32.  Hence  an  in¬ 
crease  in  the  barrier  level  of  less  than  6%  decreased  the  passage  prob¬ 
ability  by  20%. 


CHAPTER  5 


SUMMARY,  CONCLUSIONS  AND  RECOMMENDATIONS 


5.1  SUMMARY 

A  technique  to  explicitly  account  for  the  randomness  of  structural 
stiffness  and  damping  in  the  computation  of  the  structural  response  was 
developed.  The  technique  is  suitable  for  use  in  cases  of  botn  single 
and  multiple  degrees  of  freedom. 

The  technique  was  applied  to  several  problems  in  Chapters  2  and  3 
to  demonstrate  the  influence  of  these  underlying  uncertainties  on  the 
response  of  typical  structures. 

The  results  of  some  of  these  computations  were  used  to  compute  the 
probability  of  the  response  exceeding  a  given  value  using  a  new  method 
developed  in  Chapter  4. 

The  methods  of  Chapters  2  and  3  are  limited  to  linear  elastic 
structures,  however  they  may  be  readily  extended  to  cases  of  nonlinear¬ 
ity  or  plastic  response. 

The  method  for  computing  the  passage  probabilities  developed  in 
Chapter  4  may  be  applied  to  any  continuous  nonstationary  random  process 
which  satisfies  the  assumptions  made  in  Chapter  4. 

5.2  CONCLUSIONS 

The  randomness  of  a  systems  stiffness  and  damping  will  effect  the 
mean  and  standard  deviation  of  Liat  systems  response.  The  uncertainties 
will  have  only  minimal  influence  when  the  system  response  is  mean  zero. 
This  is  to  be  expected  since  the  systems  stiffness  and  damping  only 
enter  the  problem  when  the  system  is  displaced. 
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The  relative  influence  of  the  system  stiffness  appears  constant  for 
all  non-zero  deflections.  The  mean  at  any  time  may  oe  approximated  by, 

=  Xp  •  (1  +  Cov2]  (5-1) 

where  is  the  mean  of  X,  Xq  is  the  deterministic  value  and  covk 
is  the  coefficient  of  variation  for  the  stiffness. 

The  relative  influence  of  the  damping  depends  on  the  level  of  damp¬ 
ing,  and  has  no  effect  on  the  static  deflections.  For  lightly  damped 
systems,  e  <_  0.15,  the  influence  of  the  randomness  of  the  damping  may 
probably  be  ignored  for  most  problems. 

The  resuUs  of  Chapter  4  indicate  that  the  randomness  in  the  stiff¬ 
ness  must  be  accounted  for  when  a  high  degree  of  survivability  is 
required.  A  blast  excited  slab  had  a  significant  probability  of  exceed¬ 
ing  the  theoretical  maximum  vaiue  obtained  for  the  deterministic  system. 
5.3  RECOMMENDATIONS 

The  problem  of  MOOF  system  response  needs  more  work,  specifically 
an  efficient  method  of  accounting  for  the  correlation  of  the  random  var¬ 
iables  within  the  stiffness  and  damping  matrices.  This  may  be  easier  in 
a  plane  strain  or  plane  stress  computai'lon  where  the  system  is  a  con¬ 
tinuous  medium. 

A  method  for  accounting  for  sytem  uncertainties  in  cases  where  the 
response  is  not  (iiean  zero  should  d'*so  be  developed.  This  includes  suc.> 
important  cases  as  the  response  to  earthquake,  and  turbulence. 


APPENDIX  A 


PROGRAM  RSDOF 
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